Teil A: Hauspreise analysieren und einfach prognostizieren¶
Aufgabe A-0: Bibliotheken, Funktionen, Datenstrukturen¶
Alle für die weiteren Aufgaben benötigten Bibliotheken und Funktionen sind zentral an dieser Stelle einzubinden bzw. zu definieren. Ebenso sollen alle globalen Variablen in diesem Teil angegeben werden. Des Weiteren sollen Datenstrukturen für den Vergleich der Laufzeiten der Teile A, B und C sowie Datenstrukturen für den Vergleich der Laufzeiten und Modellgüten RMSE bezüglich der Validierungsstichprobe für die Modelle der Aufgaben A-7 und B-1, B-2 b) bis B-4 erstellt werden.
!pip install catboost
Requirement already satisfied: catboost in /usr/local/lib/python3.10/dist-packages (1.2.5) Requirement already satisfied: graphviz in /usr/local/lib/python3.10/dist-packages (from catboost) (0.20.3) Requirement already satisfied: matplotlib in /usr/local/lib/python3.10/dist-packages (from catboost) (3.7.1) Requirement already satisfied: numpy>=1.16.0 in /usr/local/lib/python3.10/dist-packages (from catboost) (1.25.2) Requirement already satisfied: pandas>=0.24 in /usr/local/lib/python3.10/dist-packages (from catboost) (2.0.3) Requirement already satisfied: scipy in /usr/local/lib/python3.10/dist-packages (from catboost) (1.11.4) Requirement already satisfied: plotly in /usr/local/lib/python3.10/dist-packages (from catboost) (5.15.0) Requirement already satisfied: six in /usr/local/lib/python3.10/dist-packages (from catboost) (1.16.0) Requirement already satisfied: python-dateutil>=2.8.2 in /usr/local/lib/python3.10/dist-packages (from pandas>=0.24->catboost) (2.8.2) Requirement already satisfied: pytz>=2020.1 in /usr/local/lib/python3.10/dist-packages (from pandas>=0.24->catboost) (2023.4) Requirement already satisfied: tzdata>=2022.1 in /usr/local/lib/python3.10/dist-packages (from pandas>=0.24->catboost) (2024.1) Requirement already satisfied: contourpy>=1.0.1 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (1.2.1) Requirement already satisfied: cycler>=0.10 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (0.12.1) Requirement already satisfied: fonttools>=4.22.0 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (4.51.0) Requirement already satisfied: kiwisolver>=1.0.1 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (1.4.5) Requirement already satisfied: packaging>=20.0 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (24.0) Requirement already satisfied: pillow>=6.2.0 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (9.4.0) Requirement already satisfied: pyparsing>=2.3.1 in /usr/local/lib/python3.10/dist-packages (from matplotlib->catboost) (3.1.2) Requirement already satisfied: tenacity>=6.2.0 in /usr/local/lib/python3.10/dist-packages (from plotly->catboost) (8.2.3)
from google.colab import drive
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import random
import numpy as np
import math
import sklearn.linear_model
import sklearn.preprocessing
import sklearn.metrics
import catboost
import sklearn
import time
from tensorflow import keras
import pickle
import tensorflow as tf
import os
import lightgbm as lgb
import sklearn.cluster
import xgboost as xgb
import geopandas as gpd
import psutil
seed = 43
random.seed(seed)
np.random.seed(seed)
tf.random.set_seed(seed)
# Laufzeiten als Dictionary
times = {}
# Dictionary der RMSE für Validierungsdaten und Testdaten
rmses = {}
rmsestest = {}
runtime_folder='cpu/'
drive_path='/content/drive'
folder_path=drive_path+'/MyDrive/cads/pruefung/'
results_path=folder_path+runtime_folder
def save_times(times):
timepath=results_path+'times.pkl'
with open(timepath, 'wb') as file:
pickle.dump(times, file)
def plotcali(data, column, label, cal_img) :
plot_axis = data.plot(
kind="scatter", x="longitude", y="latitude",
s=data["population"] / 100,
xlabel="Longitude", ylabel="Latitude",
label="Population",
c=column, colormap="jet", colorbar=True,
legend=True, sharex=False, figsize=(20, 10))
plot_axis.collections[0].colorbar.set_label(label)
axisrange = [-124.55, -113.95, 32.45, 42.05]
plt.axis(axisrange)
plt.imshow(cal_img, extent=axisrange)
plt.show()
def rmse(y_val, y_pred):
return np.sqrt(sklearn.metrics.mean_squared_error(y_val, y_pred))
def print_results(y_val, y_pred):
print("RMSE:",rmse(y_val, y_pred))
rmst = np.std(y_val, ddof=0)
print("RMST (zum Vergleich):",rmst)
r2 = sklearn.metrics.r2_score(y_val, y_pred)
print("R^2:",r2)
def plotRelFreq(data, column, column_label, height):
rel_freq = (data[column].value_counts() / len(data))*100.0
plotdata = pd.DataFrame({column: rel_freq.index,'relfreq': rel_freq.values,})
plt.figure(figsize=(20, height))
sns.barplot(x='relfreq', y=column, data=plotdata, orient='h')
plt.yticks(fontsize=14)
plt.xticks(fontsize=14)
plt.xlabel('Relative Häufigkeit in %', fontsize=14)
plt.ylabel(column_label, fontsize=14)
plt.tight_layout()
plt.show()
# Insgesamt log(Neuronen) Layer
# Anzahl der Neuronen pro Layer konstant
def simple_layer_shape(neurons):
layers = int(math.log(neurons))
npl = int(neurons/layers)
layer_neurons_list = [npl for i in range(0,layers)]
return layer_neurons_list
# Einfaches Regressionsnetzwerk
def simple_network(middle_neurons, num_features):
layer_shape = simple_layer_shape(middle_neurons)
layers = [keras.layers.Dense(layer_shape[0], activation='relu', input_shape=(num_features,))]
layers.extend([keras.layers.Dense(layer_size, activation='relu') for layer_size in layer_shape[1:]])
layers.append(keras.layers.Dense(1))
model = keras.Sequential(layers)
model.compile(optimizer='adam',
loss='mean_squared_error',
metrics=['mean_squared_error'])
return model
# Einfaches Regressionsnetzwerk mit Regularisierung und vorgebbarer Lernrate
def simple_regularized_network(middle_neurons, num_features, reg_strength, learning_rate=0.001):
layer_shape = simple_layer_shape(middle_neurons)
layers = [keras.layers.Dense(layer_shape[0], activation='relu', kernel_regularizer=keras.regularizers.l2(reg_strength), input_shape=(num_features,))]
layers.extend([keras.layers.Dense(layer_size, activation='relu',kernel_regularizer=keras.regularizers.l2(reg_strength)) for layer_size in layer_shape[1:]])
layers.append(keras.layers.Dense(1))
model = keras.Sequential(layers)
optimizer = keras.optimizers.Adam(learning_rate=learning_rate)
model.compile(optimizer=optimizer,
loss='mean_squared_error',
metrics=['mean_squared_error'])
return model
def save_model(model, history, path, name):
model.save(path+name+'_model.h5')
with open(path+name+'_history.pkl', 'wb') as file:
pickle.dump(history.history, file)
def load_model(path, name):
print(path+name+'_model.h5')
model = keras.models.load_model(path+name+'_model.h5')
with open(path+name+'_history.pkl', 'rb') as file:
loaded_history_dict = pickle.load(file)
history = type('History', (), {'history': loaded_history_dict})
return model, history
# Modell entweder trainieren und abspeichern oder gespeichertes Modell laden
def train_model(model, x_train, y_train, epochs, name, retrain_model=True, do_save_model=True, verbosity=0):
if retrain_model:
history = model.fit(x_train.values.astype(float), y_train.values.astype(float), epochs=epochs, batch_size=min(512,len(x_train)), validation_split=0.2, verbose=verbosity)
if do_save_model:
save_model(model, history, results_path, name)
else:
model, history = load_model(results_path, name)
return model, history
def train_model_early_stopping(model, x_train, y_train, epochs, name, retrain_model=True, do_save_model=True, verbosity=0):
if retrain_model:
early_stopping = keras.callbacks.EarlyStopping(monitor='val_mean_squared_error', patience=100, restore_best_weights=True)
history = model.fit(x_train.values.astype(float), y_train.values.astype(float), epochs=epochs, batch_size=min(512,len(x_train)), validation_split=0.2, callbacks=[early_stopping],verbose=verbosity)
if do_save_model:
save_model(model, history, results_path, name)
else:
model, history = load_model(results_path, name)
return model, history
def train_embed_model(model, x_train, y_train, epochs, name, retrain_model=True, do_save_model=True, verbosity=0):
if retrain_model:
early_stopping = keras.callbacks.EarlyStopping(monitor='val_mean_squared_error', patience=100, restore_best_weights=True)
history = model.fit([x_train['ocean_proximity'].cat.codes,x_train['county_name'].cat.codes, x_train[numeric_features].values.astype(float)], y_train.values.astype(float), epochs=epochs, batch_size=min(512,len(x_train)), validation_split=0.2, callbacks=[early_stopping],verbose=verbosity)
if do_save_model:
save_model(model, history, results_path, name)
else:
model, history = load_model(results_path, name)
return model, history
def plot_embedding(feature_name,embedding_model, x_train, title, x_size=20, y_size=10, fontsize=16):
code_to_category = {code: cat for code,cat in zip(x_train[feature_name].cat.codes, x_train[feature_name])}
embedding_predictions = embedding_model.predict(np.array(list(code_to_category.keys())))
x_coords = embedding_predictions[:,0, 0]
y_coords = embedding_predictions[:,0, 1]
labels = list(code_to_category.values())
plt.figure(figsize=(x_size, y_size))
plt.scatter(x_coords, y_coords, alpha=0.5)
for i, label in enumerate(labels):
plt.text(x_coords[i], y_coords[i], label, fontsize=fontsize, ha='left', va='bottom')
plt.xlabel('Dimension 1')
plt.ylabel('Dimension 2')
plt.title(title)
plt.show()
def plot_embed_axis(ax, x_coords, y_coords, labels,title):
ax.scatter(x_coords, y_coords, alpha=0.5)
for i, label in enumerate(labels):
ax.text(x_coords[i], y_coords[i], label, fontsize=10, ha='left', va='bottom')
ax.set_xlabel('Dim 1')
ax.set_ylabel('Dim 2')
ax.set_title(title)
def embedding_to_coords(embedding, series):
code_to_category = {code: cat for code,cat in zip(series.cat.codes, series)}
embedding_predictions = embedding.predict(np.array(list(code_to_category.keys())))
x_coords = embedding_predictions[:,0, 0]
y_coords = embedding_predictions[:,0, 1]
labels = list(code_to_category.values())
return x_coords, y_coords, labels
def ncornerpoints(polygon):
return len(polygon.exterior.coords) + sum(len(interior.coords) for interior in polygon.interiors)
retrain_models = True #Modelle entweder neu trainieren oder laden
nn_verbosity = 0
drive.mount(drive_path)
if not retrain_models:
timepath=results_path+'times.pkl'
if os.path.exists(timepath):
with open(timepath, 'rb') as file:
times = pickle.load(file)
del timepath, file
Drive already mounted at /content/drive; to attempt to forcibly remount, call drive.mount("/content/drive", force_remount=True).
Aufgabe A-1: Datensatz einlesen und prüfen¶
Der Ihnen zur Verfügung gestellte Datensatz california_housing_county_sample-de.csv enthält demografische Informationen zum US-Bundesstaat Kalifornien. Insgesamt enthält der Datensatz die folgenden 12 Merkmale: [...] Es ist zu beachten, dass es sich bei den Angaben der Merkmale 3 bis 9 um aggregierte Informationen auf Bezirksebene (Wahlbezirke 1992) handelt. Die nachfolgenden Aufgaben sind der Reihe nach durchzuführen:
Aufgabe A-1 a)¶
Der zur Verfügung gestellte Datensatz ist einzulesen. Im Anschluss sind die Anzahl an Zeilen und Spalten sowie sieben zufällige Zeilen auszugeben.
data = pd.read_csv(folder_path+'california_housing_county_sample-de.csv',delimiter=';',decimal=',',thousands='.')
rows, columns = data.shape
display('Zeilen: '+str(rows))
display('Spalten: '+str(columns))
display('7 zufällige Zeilen:')
display(data.sample(n=7, random_state=seed))
'Zeilen: 20640'
'Spalten: 12'
'7 zufällige Zeilen:'
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | ocean_proximity | county_name | sample | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 34.05 | -117.27 | 34 | 1703 | 395.0 | 849 | 359 | 3.1607 | 138200 | INLAND | San Bernardino | A |
| 4402 | 33.88 | -117.57 | 35 | 1755 | 446.0 | 1453 | 428 | 2.3160 | 119400 | INLAND | Riverside | C |
| 1929 | 33.60 | -117.90 | 25 | 2465 | 585.0 | 906 | 472 | 3.6538 | 500001 | <1H OCEAN | Orange | A |
| 11551 | 35.37 | -118.97 | 34 | 1379 | 333.0 | 1156 | 315 | 1.7197 | 48900 | INLAND | Kern | A |
| 9882 | 34.20 | -118.36 | 14 | 1878 | 614.0 | 1874 | 559 | 2.5267 | 231800 | <1H OCEAN | Los Angeles | C |
| 3896 | 33.85 | -118.39 | 24 | 4373 | 871.0 | 1830 | 824 | 5.7128 | 366200 | <1H OCEAN | Los Angeles | C |
| 16638 | 37.83 | -122.25 | 52 | 1279 | 287.0 | 534 | 291 | 3.1429 | 231600 | NEAR BAY | Alameda | B |
del drive_path
del rows, columns
Aufgabe A-1 b)¶
Die Datentypen des eingelesenen Datensatzes sind auszugeben, zu prüfen und zu kommentieren. Nominale Merkmale sollen den Datentyp für Kategorien (Python: category, R: factor) erhalten. Gegebenenfalls erforderliche Anpassungen sind durchzuführen und die finalen Datentypen sind auszugeben.
Datentypen:
data.dtypes
latitude float64 longitude float64 housing_median_age int64 total_rooms int64 total_bedrooms float64 population int64 households int64 median_income float64 median_house_value int64 ocean_proximity object county_name object sample object dtype: object
latitude, longitude, total_bedrooms und median_income sind Gleitkommazahlen.
Bei total_bedrooms sind semantisch eigentlich Ganzzahlen zu erwarten. Das ist weiter zu untersuchen.
housing_median_age, total_rooms, population, households und median_house_value sind Ganzzahlen. Falls weitere Berechnungen (Durchschnittsbildung etc.) auf diesen notwendig sind können sie zu Gleitkommazahlen konvertiert werden.
ocean_proximity, county_name und sample sind vom typ object und beinhalten Strings. Diese werden wir weiter konvertieren.
Testen auf fehlende Einträge:
for column in data.columns :
print(column," hat ",data[data[column].isna() == True].size," fehlende Einträge")
del column
latitude hat 0 fehlende Einträge longitude hat 0 fehlende Einträge housing_median_age hat 0 fehlende Einträge total_rooms hat 0 fehlende Einträge total_bedrooms hat 2484 fehlende Einträge population hat 0 fehlende Einträge households hat 0 fehlende Einträge median_income hat 0 fehlende Einträge median_house_value hat 0 fehlende Einträge ocean_proximity hat 0 fehlende Einträge county_name hat 0 fehlende Einträge sample hat 0 fehlende Einträge
total_bedrooms hat fehlende Einträge. Das könnte der Grund für die Gleitkommakonvertierung sein (int wäre als 0 nicht eindeutig erkennbar, float als NaN aber schon). Die Bereinigung machen wir später.
Die kategoriellen Daten werden in den Datentyp "category" konvertiert.
categorical_columns = ['ocean_proximity','county_name','sample']
for merkmal in categorical_columns:
data[merkmal] = data[merkmal].astype('category')
del merkmal
Nach den Anpassungen haben wir folgende Datentypen:
data.dtypes
latitude float64 longitude float64 housing_median_age int64 total_rooms int64 total_bedrooms float64 population int64 households int64 median_income float64 median_house_value int64 ocean_proximity category county_name category sample category dtype: object
Aufgabe A-1 c)¶
Das Merkmal sample, das fest vorgegebene Stichproben für das spätere Trainieren (A), Validieren (B) und Testen (C) der Modelle definiert, ist auszuzählen und die Stichprobengrößen sind zu kommentieren.
Anzahl der Stichproben pro Sample:
data['sample'].value_counts()
sample A 14421 B 3110 C 3109 Name: count, dtype: int64
Anteile an der Gesamtanzahl der Stichproben
data['sample'].value_counts() / data['sample'].size
sample A 0.698692 B 0.150678 C 0.150630 Name: count, dtype: float64
A umfasst mit fast 70% die meisten Stichproben, B und C haben fast eine ähnliche Anzahl mit jeweils knapp über 15%. A eigent sich damit als Trainingssample, während B und C als Validierungs- bzw. Testsample verwendet werden können (insofern sich im weiteren Verlauf nicht herausstellt, dass die Samples einen Bias besitzen, also nicht zufällig zugeteilt wurden)
Aufgabe A-2: Numerische Merkmale analysieren, aufbereiten und visualisieren¶
Aufgabe A-2 a)¶
Für die numerische Merkmale sind mindestens die Kennzahlen count, mean,
std, min/max sowie die das 25%-, 50%- und 75%-Quantil auszugeben. Zudem soll ermittelt werden, bei welchen Merkmalen fehlende Werte (und falls ja, in welcher absoluten und prozentualen Größenordnung) vorliegen.
Spalten mit Zahlenwerten:
numeric_columns = [column for column, dtype in data.dtypes.items() if dtype in ['int64', 'float64']]
display(numeric_columns)
['latitude', 'longitude', 'housing_median_age', 'total_rooms', 'total_bedrooms', 'population', 'households', 'median_income', 'median_house_value']
Übersicht der Zahlenwerte
data[numeric_columns].describe()
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | |
|---|---|---|---|---|---|---|---|---|---|
| count | 20640.000000 | 20640.000000 | 20640.000000 | 20640.000000 | 20433.000000 | 20640.000000 | 20640.000000 | 20640.000000 | 20640.000000 |
| mean | 35.631861 | -119.569704 | 28.639486 | 2635.763081 | 537.870553 | 1425.476744 | 499.539680 | 3.870671 | 206855.816909 |
| std | 2.135952 | 2.003532 | 12.585558 | 2181.615252 | 421.385070 | 1132.462122 | 382.329753 | 1.899822 | 115395.615874 |
| min | 32.540000 | -124.350000 | 1.000000 | 2.000000 | 1.000000 | 3.000000 | 1.000000 | 0.499900 | 14999.000000 |
| 25% | 33.930000 | -121.800000 | 18.000000 | 1447.750000 | 296.000000 | 787.000000 | 280.000000 | 2.563400 | 119600.000000 |
| 50% | 34.260000 | -118.490000 | 29.000000 | 2127.000000 | 435.000000 | 1166.000000 | 409.000000 | 3.534800 | 179700.000000 |
| 75% | 37.710000 | -118.010000 | 37.000000 | 3148.000000 | 647.000000 | 1725.000000 | 605.000000 | 4.743250 | 264725.000000 |
| max | 41.950000 | -114.310000 | 52.000000 | 39320.000000 | 6445.000000 | 35682.000000 | 6082.000000 | 15.000100 | 500001.000000 |
Auffällig ist, dass median_incom offensichtlich nicht direkt über Doller-Werte repräsentiert ist. Möglicherweise sind dies Mehrfache von 10.000$. Das glatte Maximum deutet darauf hin, dass hier ab einer Obergrenze abgeschnitten bzw. minimiert wurde.
Fehlende Werte:
for column in data.columns :
nas = data[data[column].isna() == True].size
if(nas>0):
print(column," hat ",nas," fehlende Einträge")
print("Das ist ein Anteil von ",nas / data[column].size," des Datensatzes.")
del column
del nas
total_bedrooms hat 2484 fehlende Einträge Das ist ein Anteil von 0.12034883720930233 des Datensatzes.
Es fehlen also etwas mehr als 12% der Werte für total_bedrooms.
Aufgabe A-2 b)¶
Anhand des Datensatzes ist zu begründen und geeignet grafisch zu visualisieren, warum eine Ersetzung der fehlenden Werte bei total_bedrooms durch den Median nicht empfehlenswert ist (mindestens zwei konkrete Zeilen des Datensatzes sind als Beispiele auszugeben) und was für einen Bezug auf total_rooms spricht. Im Anschluss sind für das Feature total_bedrooms die fehlenden Werte durch den mittleren Quotienten zwischen total_bedrooms und total_rooms (berechnet auf Basis der Datensätze ohne fehlende Werte) zu ersetzen. Die Ersetzung ist geeignet zu verifizieren.
sns.histplot(data=data, x='total_bedrooms', bins=10)
plt.title("Historgramm")
plt.show()
Im Histogram ist zu sehen, dass total_bedrooms eine sehr schiefe Verteilung besitzt. Das spricht zunächst gegen den Median als Ersetzungswert.
sns.boxplot(data=data, x='total_bedrooms')
plt.title("Boxplot")
plt.show()
display(data['total_bedrooms'].describe())
count 20433.000000 mean 537.870553 std 421.385070 min 1.000000 25% 296.000000 50% 435.000000 75% 647.000000 max 6445.000000 Name: total_bedrooms, dtype: float64
Der Boxplot bestätigt die schiefe Verteilung. Mit einem Minimum von 1 und einem Maximum von 6445 gibt es gegenüber den Quantilen auch extreme Werte in den Tails, insbesondere im oberen Quartil.
data['total_bedrooms'].median()
435.0
print("Minimum Population:")
display(data.loc[data['population'].idxmin()])
Minimum Population:
latitude 34.04 longitude -118.44 housing_median_age 16 total_rooms 18 total_bedrooms 6.0 population 3 households 4 median_income 0.536 median_house_value 350000 ocean_proximity <1H OCEAN county_name Los Angeles sample B Name: 6895, dtype: object
print("Maximum Population:")
display(data.loc[data['population'].idxmax()])
Maximum Population:
latitude 33.35 longitude -117.42 housing_median_age 14 total_rooms 25135 total_bedrooms 4819.0 population 35682 households 4769 median_income 2.5729 median_house_value 134400 ocean_proximity <1H OCEAN county_name San Diego sample C Name: 1705, dtype: object
Eine Nutzung des Medians von 435 wäre insbesondere in den Extremen der Population komplett unplausibel.
sns.scatterplot(data=data, x='total_bedrooms', y='total_rooms')
plt.show()
data['total_bedrooms'].corr(data['total_rooms'])
0.9303795046865072
Es besteht offensichtlich eine starke Korrelation zwischen Schlafzimmern und Räumen, sie sind nahezu linear zueinander. Das zeigt auch der Korrelationskoeffizient.
bedrooms_per_room = (data[data['total_bedrooms'].isna() == False]['total_bedrooms'] / data[data['total_bedrooms'].isna()==False]['total_rooms']).mean()
print("Mittlerer Anteil Schlafzimmer an allen Räumen: ", bedrooms_per_room)
Mittlerer Anteil Schlafzimmer an allen Räumen: 0.2130388304808513
Fehlende Werte der Schlafzimmer werden via Faktor aus den Räumen ersetzt
missing_indices_sample = data[data['total_bedrooms'].isna()].sample(n=5, random_state=seed).index
print("Sample der Daten vor Imputation")
display(data.loc[missing_indices_sample])
Sample der Daten vor Imputation
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | ocean_proximity | county_name | sample | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 6446 | 34.01 | -117.86 | 16 | 4632 | NaN | 3038 | 727 | 5.1762 | 264400 | <1H OCEAN | Los Angeles | A |
| 4212 | 33.87 | -117.92 | 33 | 1597 | NaN | 1888 | 423 | 3.0550 | 157800 | <1H OCEAN | Orange | B |
| 7594 | 34.07 | -118.37 | 50 | 2519 | NaN | 1117 | 516 | 4.3667 | 405600 | <1H OCEAN | Los Angeles | A |
| 7407 | 34.06 | -118.28 | 42 | 2472 | NaN | 3795 | 1179 | 1.2254 | 162500 | <1H OCEAN | Los Angeles | A |
| 17686 | 38.01 | -120.37 | 30 | 473 | NaN | 242 | 93 | 2.5417 | 123200 | INLAND | Tuolumne | C |
Fehlende Werte der Schlafzimmer werden via Faktor aus den Räumen ersetzt
data['total_bedrooms'] = data['total_bedrooms'].fillna(data['total_rooms']*bedrooms_per_room)
del bedrooms_per_room
print("Sample nach Imputation")
display(data.loc[missing_indices_sample])
del missing_indices_sample
Sample nach Imputation
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | ocean_proximity | county_name | sample | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 6446 | 34.01 | -117.86 | 16 | 4632 | 986.795863 | 3038 | 727 | 5.1762 | 264400 | <1H OCEAN | Los Angeles | A |
| 4212 | 33.87 | -117.92 | 33 | 1597 | 340.223012 | 1888 | 423 | 3.0550 | 157800 | <1H OCEAN | Orange | B |
| 7594 | 34.07 | -118.37 | 50 | 2519 | 536.644814 | 1117 | 516 | 4.3667 | 405600 | <1H OCEAN | Los Angeles | A |
| 7407 | 34.06 | -118.28 | 42 | 2472 | 526.631989 | 3795 | 1179 | 1.2254 | 162500 | <1H OCEAN | Los Angeles | A |
| 17686 | 38.01 | -120.37 | 30 | 473 | 100.767367 | 242 | 93 | 2.5417 | 123200 | INLAND | Tuolumne | C |
print("Noch fehlende Werte für Schlafzimmer: ",data[data['total_bedrooms'].isna()].size)
Noch fehlende Werte für Schlafzimmer: 0
Aufgabe A-2 c)¶
Zu den bereits vorhandenen Merkmalen sollen noch die folgenden Merkmale zum Datensatz hinzugefügt werden: • bedrooms_per_house = total_bedrooms / households • rooms_per_house = total_rooms / households Die Berechnung ist an zwei ausgewählten Datensätzen zu demonstrieren.
Ergänzung von Schlafzimmer pro Haushalt und Räume pro Haushalt.
data['bedrooms_per_house'] = data['total_bedrooms'] / data['households']
data['rooms_per_house'] = data['total_rooms'] / data['households']
if 'bedrooms_per_house' not in numeric_columns: numeric_columns.append('bedrooms_per_house')
if 'rooms_per_house' not in numeric_columns: numeric_columns.append('rooms_per_house')
data.sample(n=2, random_state=seed)
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | ocean_proximity | county_name | sample | bedrooms_per_house | rooms_per_house | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 34.05 | -117.27 | 34 | 1703 | 395.0 | 849 | 359 | 3.1607 | 138200 | INLAND | San Bernardino | A | 1.100279 | 4.743733 |
| 4402 | 33.88 | -117.57 | 35 | 1755 | 446.0 | 1453 | 428 | 2.3160 | 119400 | INLAND | Riverside | C | 1.042056 | 4.100467 |
Aufgabe A-2 d)¶
Die numerischen Merkmale sind zusammen mit einem Kerndichteschätzer zu visualisieren und die Ausgaben zu kommentieren.
Kerndichteplot für jedes numerisches Merkmal
for column in numeric_columns:
plt.figure(figsize=(20, 4))
sns.histplot(data=data[column], bins='auto', kde=True, fill=True)
plt.show()
del column
Es fällt auf, dass die Längen- und Breitengrade zwei Häufungen besitzen.
Bei den Attributen, die mit der Population skalieren (population, total_rooms, total_bedrooms, households) hat man wieder eine stark asymetrische Verteilung mit einem dünnen Peak bei niedrigen Werten und einem sehr langen Tail bei größeren Werten. Einkommen und Hauswert haben zwar eine ähnliche Tendenz, sind aber nicht so extrem asymetrisch wie die Populationsgrößen. Das Einkommen ist wie weiter oben beschrieben nicht in USD angegeben, sondern vorskaliert. Das mittlere Alter der Häuser ist mehr oder weniger symmetrisch um die 25 Jahre verteilt.
Die Daten des Median-Hauswertes scheinen mit 50001$ minimiert worden zu sein. Ein ähnliches "Capping" wurde scheinbar auch beim Median-Einkommen und beim Median-Hausalter gemacht. Je nach weiterem Anwendungsfall kann es notwendig sein, diese Daten zu bereinigen.
bedrooms_per_house und rooms_per_house haben ebenfalls einen starken Peak bei niedrigen Werten und einen langen Tail in den oberen Bereichen. Hier ist denkbar dass es lediglich ein paar wenige Ausreißer sind, welche die größeren Werte beschreiben (Gebiete mit fast ausschließlich Hotels o.ä.)
Aufgabe A-2 e)¶
Die Korrelationen der numerischen Merkmale sind zu berechnen und zu visualisieren. Im Anschluss ist auf Basis dieser Ergebnisse zu begründen, warum die Merkmale households, total_rooms, total_bedrooms und bedrooms_per_house aus dem Datensatz entfernt werden können. Das ist im Anschluss durchzuführen. Im Anschluss sind die Korrelationen der verbleibenden Merkmale erneut zu plotten.
Darstellung der Korrelationen
data[numeric_columns].corr()
| latitude | longitude | housing_median_age | total_rooms | total_bedrooms | population | households | median_income | median_house_value | bedrooms_per_house | rooms_per_house | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| latitude | 1.000000 | -0.924664 | 0.011173 | -0.036100 | -0.066127 | -0.108785 | -0.071035 | -0.079809 | -0.144160 | 0.070192 | 0.106389 |
| longitude | -0.924664 | 1.000000 | -0.108197 | 0.044568 | 0.068416 | 0.099773 | 0.055310 | -0.015176 | -0.045967 | 0.012916 | -0.027540 |
| housing_median_age | 0.011173 | -0.108197 | 1.000000 | -0.361262 | -0.321501 | -0.296244 | -0.302916 | -0.119034 | 0.105623 | -0.078659 | -0.153277 |
| total_rooms | -0.036100 | 0.044568 | -0.361262 | 1.000000 | 0.930847 | 0.857126 | 0.918484 | 0.198050 | 0.134153 | 0.030251 | 0.133798 |
| total_bedrooms | -0.066127 | 0.068416 | -0.321501 | 0.930847 | 1.000000 | 0.877492 | 0.978765 | -0.005590 | 0.051377 | 0.046766 | 0.002541 |
| population | -0.108785 | 0.099773 | -0.296244 | 0.857126 | 0.877492 | 1.000000 | 0.907222 | 0.004834 | -0.024650 | -0.066860 | -0.072213 |
| households | -0.071035 | 0.055310 | -0.302916 | 0.918484 | 0.978765 | 0.907222 | 1.000000 | 0.013033 | 0.065843 | -0.055812 | -0.080598 |
| median_income | -0.079809 | -0.015176 | -0.119034 | 0.198050 | -0.005590 | 0.004834 | 0.013033 | 1.000000 | 0.688075 | -0.057669 | 0.326895 |
| median_house_value | -0.144160 | -0.045967 | 0.105623 | 0.134153 | 0.051377 | -0.024650 | 0.065843 | 0.688075 | 1.000000 | -0.044737 | 0.151948 |
| bedrooms_per_house | 0.070192 | 0.012916 | -0.078659 | 0.030251 | 0.046766 | -0.066860 | -0.055812 | -0.057669 | -0.044737 | 1.000000 | 0.849030 |
| rooms_per_house | 0.106389 | -0.027540 | -0.153277 | 0.133798 | 0.002541 | -0.072213 | -0.080598 | 0.326895 | 0.151948 | 0.849030 | 1.000000 |
Korrelationen via Heatmap
sns.heatmap(data[numeric_columns].corr(), annot=True, fmt=".2f", cmap='coolwarm')
plt.show()
Es ist auffällig, dass die Werte total_rooms, total_bedrooms, households und population fast komplett linear zueinander sind. Aus diesem Grund können households, total_rooms und total_bedrooms aus dem Datensatz entfernt werden, da ihre Aussagekraft fast komplett durch population abgedeckt ist.
Ebenso sind bedrooms_per_house und rooms_per_house dermaßen stark korreliert dass auf bedrooms_per_house verzichtet werden kann, da dessen Aussagekraft durch rooms_per_house abgedeckt ist.
Entfernen von households, total_rooms, total_bedrooms und bedrooms_per_house:
columns_to_drop = [column for column in ['households', 'total_rooms', 'total_bedrooms','bedrooms_per_house'] if column in data ]
data.drop(columns=columns_to_drop, inplace=True)
for column in columns_to_drop:
numeric_columns.remove(column)
del column
del columns_to_drop
Korrelationen mit reduzierten Daten
sns.heatmap(data[numeric_columns].corr(), annot=True, fmt=".2f", cmap='coolwarm')
plt.show()
Aufgabe A-3: Hauspreise sowie regionale Strukturen untersuchen und visualisieren¶
Aufgabe A-3 a)¶
Visualisieren Sie den Datensatz mit Hilfe der Karte von Kalifornien (siehe: https://github.com/ageron/handson-ml3/blob/main/02_end_to_end_machine_learning_project.ipynb) und den folgenden Merkmalen: • median_house_value • housing_median_age • median_income Die Ergebnisse sind zu kommentieren und zu interpretieren.
Bild von Californien und Code für den Plot auf der Karte aus dem Beispielnotebook in der Aufgabenstellung übernommen bzw. daraus angepasst.
california_img_path = folder_path+'california.png'
cal_img = plt.imread(california_img_path)
Plot für den Median Wert eines Hauses (median_house_value)
plotcali(data, column="median_house_value", label="Median House Value (USD)", cal_img=cal_img)
Man sieht direkt die Aggregation an hohen Hauspreisen an den Ballungszentren rund um
- San Franscisco, San Jose und der generellen Bay-Area (ca. 38 Latitude, -122 Longitude)
- Los Angeles (ca. 34 Lat, -118 Long) und dem größeren Umkreis an der Küste mit Santa Barbara und San Diego
Generell besitzen hauptsächlich die Küstengebiete hohe Hauspreise. Die Hauspreise im Inland sind deutlich niedriger. Ausnahmen bilden große Städte wie Sacramento und dem Gebiet rund um Lake Tahoe an der Grenze zu Nevada.
Plot für das Medianalter eines Hauses (housing_median_age)
plotcali(data, column="housing_median_age", label="Median House Age", cal_img=cal_img)
Die Medianalter der Häuser sind (offensichtlich) gleichmäßiger verteilt als die Hauswerte. Interessant sind vor allem die Gebiete mit hoher Population und niedrigem Hausalter, die auf eine hohe Rate von Neubauten schließen lassen. Das ist insbesondere an einem Gürtel rund um Los Angeles und in Sacramento zu beobachen. Die schon lange bestehenden großen Stadtzentren haben naturgemäß auch ein hohes Hausalter.
Plot für das Medianeinkommen (median_income)
plotcali(data, column="median_income", label="Median Income (Scaled)",cal_img=cal_img)
Wie man erwarten könnte decken sich Gebiete mit hohem Medianeinkommen und Gebiete mit hohen Hauspreisen.
del cal_img
del california_img_path
Aufgabe A-3 b)¶
Betrachtet werden soll die Abhängigkeit des Hauspreises von den verfügbaren numerischen Merkmalen. Dazu ist eine geeignete Visualisierungsmöglichkeit zu wählen, auszuführen und zu interpretieren.
Mit einer Heatmap wurden in den vorherigen Aufgabenteilen bereits die Korrelationen ermittelt. Via Scatterplots werden nun alle numerischen Merkmale dem Median-Hauswert gegenübergestellt.
for column in [column for column in numeric_columns if column != 'median_house_value']:
sns.scatterplot(x=column, y='median_house_value', data=data)
plt.show()
del column
Latitude und Longitude zeigen hier klare Spitzen -- dies ist aber nicht überraschend, da wir den Zusammenhang zwischen Hauspreisen und den Ballungszentren der Küste bereits in der vorherigen Aufgabe gesehen haben. Das Hausalter scheint ziemlich gleichverteilt zu sein und besitzt alleine nur eine geringe Aussagekraft über den Hauswert. Aus dem Zusammenhang mit der Population lässt sich auch keine lineare Abhängigkeit herauslesen. Das Medianeinkommen scheint der stärkste Prädiktor für den Medianhauswert zu sein, hier gibt es eine klare Korrelation. Die Räume pro Haus scheinen wenn dann nur einen sehr schwachen linearen Zusammenhang mit den Hauswerten zu besitzen.
Aufgabe A-4: Erste Hauspreisprognosemodelle auf Basis der numerischen Features¶
Ziel dieses Abschnitts ist die Erstellung erster Modelle basierend rein auf den numerischen Merkmalen.
Aufgabe A-4 a)¶
Es ist anzugeben, wie die Fehler- bzw. Gütemaße Root Mean Squared Error (RMSE) bzw. R² eines Prognoseverfahrens bei einem Regressionsproblem berechnet werden. Im Anschluss soll der mathematische Zusammenhang zwischen den beiden Größen aufgezeigt werden. Dabei soll einerseits darauf eingegangen werden, welche Rolle ein naives Nullmodell einnimmt, das stets den Mittelwert des Zielmerkmals ausgibt, und andererseits begründet werden, wieso der RMSE bei steigendem R² sinkt und bei fallendem R² steigt
In der üblichen Notation für Regressionsprobleme mit Zielgröße $y$, den tatsächlichen Zielwerten $(y_1,...,y_n)$ und den vorhergesagten Werten $(\hat{y}_1,...,\hat{y}_n)$ sind $$\text{RMSE} = \sqrt{\frac{1}{n}\sum_{i=1}^{n}(y_i - \hat{y}_i)^2}$$ und $$R^2 = 1 - \frac{\sum_{i=1}^{n}(y_i - \hat{y}_i)^2}{\sum_{i=1}^{n}(y_i - \bar{y})^2}$$ wobei $\bar{y}$ der Mittelwert der tatsächlichen Zielgrößen ist.
Mit $$\text{RMST}=\sqrt{\frac{1}{n}\sum_{i=1}^{n}(y_i - \bar{y})^2}$$ kann man $R^2$ auch darstellen als $$R^2 = 1 - \left(\frac{\text{RMSE}}{\text{RMST}}\right)^2$$
RMSE stellt ein Fehlermaß dar und bestraft die Abweichungen der Vorhersagen in absoluter Höhe. Es berücksichtigt die vorhandene Varianz der Zielgrößen nicht, sodass es für die Bewertung des Fehlers keinen Unterschied macht, ob der Datensatz eine hohe Varianz aufweist. EIn naives Nullmodell mit Schätzwert $\bar{y}$ hätte damit einen Fehlerwert von $$\text{RMST}=\sqrt{\frac{1}{n}\sum_{i=1}^{n}(y_i - \bar{y})^2}$$ erhalten. Der Fehlerwert des Nullmodells steigt also mit der Varianz der Stichprobe.
$R^2$ hingegen ist ein Gütemaß und berücksichtigt nur das Verhältnis zwischen RMSE und RMST. Es sagt also aus, wie gut der RMSE im Verhältnis zum RMST minimiert werden kann. Ds zeigt sich auch beim naiven Nullmodell mit Schätzwert $\bar{y}$, denn dann wäre $\text{RMSE}=\text{RMST}$ und somit $R^2=0$. Die Vorhersage kann also gegenüber dem RMST nicht verbessert werden. Anders gesagt: $R^2$ gibt an, wie sehr die Vorhersage unter Nutzung der erklärenden Variablen gegenüber dem trivialen Modell verbessert werden kann.
Wegen $\text{RMSE,RMST}\geq0$ und weil RMST für eine gegebene Stichprobe konstant ist folgt aus der Darstellung $$R^2 + \left(\frac{\text{RMSE}}{\text{RMST}}\right)^2 = 1$$ trivialerweise, dass $R^2$ streng monoton fallend in $RMSE$ ist und umgekehrt $RMSE$ bei steigendem $R^2$ fallen muss.
Aufgabe A-4 b)¶
Der Datensatz ist in einen Trainingsteil (konkret: x_train & y_train) und einen Validierungsteil (konkret: x_val & y_val) zu splitten. Dabei ist die Spalte sample zu berücksichtigen: Zeilen mit dem Eintrag „A“ sollten für das Training und Zeilen mit dem Eintrag „B“ für die Validierung verwendet werden. Die entsprechend benötigten Datenstrukturen sind anzulegen.
Extraktion der Trainings- und Validierungsdaten
datanum = data[numeric_columns]
datanum_train = datanum[data['sample'] == 'A']
x_train = datanum_train.drop(['median_house_value'], axis=1)
y_train = datanum_train['median_house_value']
del datanum_train
datanum_val = datanum[data['sample'] == 'B']
x_val = datanum_val.drop(['median_house_value'], axis=1)
y_val = datanum_val['median_house_value']
del datanum_val
del datanum
Aufgabe A-4 c)¶
Für die Ergebnisse dieser und der nachfolgenden Modelle in Teil A ist eine geeignete Datenstruktur anzulegen. Diese soll für jedes Modell den Wert für den Root Mean Squared Error (RMSE) enthalten.
Die RMSEs pro Modell werden in einem Dictionary gespeichert.
display(rmses)
type(rmses)
{}
dict
Aufgabe A-4 d)¶
Zu verwenden ist ein einfaches lineares Modell ohne besondere Anpassungen. Ein Ergebnisbericht ist auszugeben und die Ergebniswerte sind zu speichern und zu kommentieren.
Fit einer linearen Regression mit Default-Einstellungen
model_LR = sklearn.linear_model.LinearRegression()
model_LR.fit(x_train, y_train)
y_pred_LR = model_LR.predict(x_val)
Ergebnisse:
print_results(y_val, y_pred_LR)
rmses['Numeric Linear'] = rmse(y_val, y_pred_LR)
del y_pred_LR
RMSE: 74508.93978488311 RMST (zum Vergleich): 114008.14636807823 R^2: 0.5728848285783494
print("Coefficients:")
print("Intercept:",model_LR.intercept_)
for i, coeff in enumerate(model_LR.coef_):
print("Feature",x_train.columns[i],":",coeff)
del i
del coeff
Coefficients: Intercept: -3885263.9207639927 Feature latitude : -45738.09430260581 Feature longitude : -46388.941936147545 Feature housing_median_age : 916.6476488154744 Feature population : -1.0334605170355644 Feature median_income : 36866.79377940718 Feature rooms_per_house : 1322.6461721243422
Das lineare Modell auf numerischen Werten hat zwar eine gewisse Vorhersagekraft, hat aber mit einem $R^2$ von unter 0.6 trotzdem noch viel Verbesserungspotenzial.
del model_LR
Aufgabe A-4 e)¶
Analog zu Teilaufgabe d) ist CatBoost mit Standardparametern zu verwenden, wobei die Ausgabe der einzelnen Iterationen sowohl hier als auch in nachfolgenden Aufgaben zu unterdrücken ist. Ein Ergebnisbericht samt allen verwendeten Parametern ist auszugeben und die Ergebniswerte sind zu speichern. Zusätzlich ist die Feature Importance auszugeben und zusammen mit den Ergebniswerten zu kommentieren.
Fit eines CatBoost mit Default-Einstellungen
model_CB = catboost.CatBoostRegressor(verbose=False)
model_CB.fit(x_train, y_train)
y_pred_CB = model_CB.predict(x_val)
Ergebnisse:
print_results(y_val, y_pred_CB)
rmses['Numeric Catboost'] = rmse(y_val, y_pred_CB)
del y_pred_CB
RMSE: 48912.79736396997 RMST (zum Vergleich): 114008.14636807823 R^2: 0.8159341447366703
Feature importance:
print("Feature Importances:")
for i, importance in enumerate(model_CB.feature_importances_):
print("Feature",x_train.columns[i],":",importance)
del i, importance
Feature Importances: Feature latitude : 29.144908189740867 Feature longitude : 28.9540488643117 Feature housing_median_age : 5.8741119369703245 Feature population : 1.7413334210991498 Feature median_income : 29.724008478316748 Feature rooms_per_house : 4.561589109561224
Hyperparameter:
model_CB.get_params()
{'loss_function': 'RMSE', 'verbose': False}
Verwendete Features
model_CB.feature_names_
['latitude', 'longitude', 'housing_median_age', 'population', 'median_income', 'rooms_per_house']
CatBoost konnte gegenüber der linearen Regression ein deutlich besseres Ergebnis erzeilen. RMSE ist deutlich gesunken und $R^2$ gestiegen. Anhand der Feature Importance sieht man, dass vor allem der Ort und das Median-Einkommen besonders in Betracht gezogen werden, was unseren bisherigen Beobachtungen und Vermutungen entspricht.
del model_CB
del x_train
del x_val
del y_train
del y_val
Aufgabe A5: Kategorielle Merkmale analysieren und aufbereiten¶
Hinweis: Der Datensatz ist für diesen Abschnitt (hinsichtlich der Features) auf die kategoriellen Merkmale zu beschränken.
Aufgabe A-5 a)¶
Zu erstellen ist eine grafische prozentuale Häufigkeitsverteilung für die nominalen Merkmale ocean_proximity und county_name. Die Ergebnisse sind zu erläutern.
Darstellung der relativen Häufigkeiten der kategoriellen Daten
plotRelFreq(data,'ocean_proximity', 'Entfernung zum Meer', 10)
Es gibt kaum Daten zu Häusern auf Inseln. Über 65% der Daten betreffen Häuser mit weniger als einer Stunde zum Meer, davon ca. 11% in der Nähe von Buchten und ca. 13% direkt in der Nähe des Ozeans.
plotRelFreq(data,'county_name', 'Landkreis', 16)
print("Los Angeles Population:",data[data['county_name']=='Los Angeles']['population'].sum())
print("San Francisco Population:", data[data['county_name']=='San Francisco']['population'].sum())
Los Angeles Population: 8779165 San Francisco Population: 702492
Auffällig ist wie häufig Los Angeles im Datensatz vertreten ist (ca. 28%). Das ist bei der hohen Population allerdings nicht überraschend. Die große Metropole San Francisco ist im Vergelich weit weniger vertreten. Der Grund dafür ist jedoch, dass sich die Bevölkerung der Region um San Francisco weiterflächiger auf mehrere Landkreise erstreckt, dazu lässt sich z,B. Santa Clara zählen, welches selbst mehr Datensätze als San Francisco besitzt.
Aufgabe A-5 b)¶
Der Einfluss der nominalen Merkmale auf die Hauspreise ist mit geeigneten Mitteln zu visualisieren. Auffälligkeiten sind zu kommentieren.
Darstellung der Auswirkungen der kategoriellen Merkmale auf den Median-Hauspreis mittels Violinplot
Hier: Violinplot in Abhängigkeit zu Ozeannähe
plt.figure(figsize=(20, 15))
sns.violinplot(x='ocean_proximity', y='median_house_value', data=data)
plt.xlabel('Nähe zum Ozean')
plt.ylabel('Median-Hauswert')
plt.grid(axis='y')
plt.show()
Häuser auf Inseln haben einen durchschnittlich sehr hohen Wert. Häuser in der Nähe der Bucht und dem Ozean haben eine breite Spanne an Werten, mit Tendenz zu hohem Wert. Bei Häusern im Inland ist die Varianz geringer, mit den niedrigsten Werten der Kategorien.
Violinplot in Abhängigkeit zum Landkreis
plt.figure(figsize=(15, 80))
sns.violinplot(x='median_house_value', y='county_name', data=data,orient='h',split=True)
plt.ylabel('Landkreis', fontsize=14)
plt.xlabel('Median-Hauswert', fontsize=14)
plt.yticks(fontsize=12)
plt.tight_layout()
plt.grid(axis='x')
plt.show()
Die Varianz der Median-Hauswerte ist von Landkreis zu Landkreis sehr unterschiedlich. Santa Barbara z.B. deckt eine sehr breite Spanne von Hauspreisen ab -- hier lässt sich mit der Information über den Landkreis nur wenig aussagen. Auf der anderen Seite hat Modoc eine ziemlich kompakte Verteilung, bei der alleine das Wissen über den Landkreis eine ziemlich gute Vorhersage über den Median-Hauswert ermöglicht.
Aufgabe A-6: Hauspreisprognosemodelle auf Basis der kategoriellen Features erstellen¶
Ziel dieser Aufgabe ist die Erstellung weiterer Modelle basierend rein auf den kategoriellen Merkmalen.
Aufgabe A-6 a)¶
Die Datenstrukturen für Training und Validierung sind analog zu Aufgabe A-4 b) anzulegen und zu verwenden.
Aufteilung in Trainingsdaten und Validierungsdaten für kategorielle Features analog zur Aufteilung bei numerischen Features. Da lineare Regression ein Encoding erfordert, wird ein One-Hot-Encoding angewendet.
categorical_features = categorical_columns.copy()
categorical_features.remove('sample')
encoded_data=pd.get_dummies(data, columns=categorical_features, drop_first=True)
columns_removed=data.columns.difference(encoded_data.columns)
columns_added=encoded_data.columns.difference(data.columns)
data_train = encoded_data[encoded_data['sample'] == 'A']
x_train = data_train[columns_added]
y_train = data_train['median_house_value']
data_val = encoded_data[encoded_data['sample'] == 'B']
x_val = data_val[columns_added]
y_val = data_val['median_house_value']
del data_train
del columns_removed
del columns_added
del data_val
del encoded_data
Aufgabe A-6 b)¶
Zu verwenden ist ein einfaches lineares Modell ohne besondere Anpassungen. Ein Ergebnisbericht ist auszugeben und die Ergebniswerte sind zu speichern und zu kommentieren. Was ergibt sich im Vergleich mit den Modellen aus Aufgabe A-4?
Lineares Regressionsmodell mit nur kategoriellen Daten
model_LR = sklearn.linear_model.LinearRegression()
model_LR.fit(x_train, y_train)
y_pred_LR = model_LR.predict(x_val)
print_results(y_val, y_pred_LR)
rmses['Categorical Linear'] = rmse(y_val, y_pred_LR)
del y_pred_LR
del model_LR
RMSE: 91193.2953907611 RMST (zum Vergleich): 114008.14636807823 R^2: 0.36018554106223566
Das lineare Modell mit nur kategoriellen Features schneidet mit einem RSME von ca. 91000 deutlich schlechter ab als das Modell mit numerischen Features und ist damit nur etwas besser als der RMST. Gleichzeitig hat es trotzdem noch eine gewisse Aussagekraft angesichts der begrenzten Information die dem Modell zur Verfügung steht. Dies liegt vermutlich, wie oben gesehen. an jenen Landkreisen welche nur eine sehr kleine Spanne von Median-Hauswerten umfassen (z.B. Modoc).
del x_train
del x_val
del y_train
del y_val
Aufgabe A-6 c)¶
Analog zu Aufgabe A-4 ist CatBoost ohne weitere Anpassungen zu verwenden. Ein Ergebnisbericht ist auszugeben und die Ergebniswerte sind zu speichern. Zusätzlich ist die (CatBoost-interne) Feature Importance auszugeben und zusammen mit den Ergebniswerten zu kommentieren. Was ergibt sich im Vergleich zu dem Modell aus Aufgabe A-4?
Hier wird absichtlich für CatBoost ein Datensatz ohne Encoding verwendet, weil das manuelle Encoding in Aufgabenteil B geschehen soll.
data_train = data[data['sample'] == 'A']
x_train = data_train[categorical_features]
y_train = data_train['median_house_value']
del data_train
data_val = data[data['sample'] == 'B']
x_val = data_val[categorical_features]
y_val = data_val['median_house_value']
del data_val
Catboost-Modell mit nur kategoriellen Daten
model_CB = catboost.CatBoostRegressor(verbose=False,cat_features=categorical_features,random_state=seed)
model_CB.fit(x_train, y_train)
y_pred_CB = model_CB.predict(x_val)
print_results(y_val, y_pred_CB)
rmses['Categorical Catboost'] = rmse(y_val, y_pred_CB)
del y_pred_CB
display(model_CB.get_params())
model_CB.feature_names_
print("Feature Importances:")
for i, importance in enumerate(model_CB.feature_importances_):
print("Feature",x_train.columns[i],":",importance)
del i, importance
del model_CB
RMSE: 89735.68699792186 RMST (zum Vergleich): 114008.14636807823 R^2: 0.38047531608783414
{'loss_function': 'RMSE',
'verbose': False,
'random_state': 43,
'cat_features': ['ocean_proximity', 'county_name']}
Feature Importances: Feature ocean_proximity : 48.94639382344655 Feature county_name : 51.05360617655345
Hier erreicht CatBoost einen nur unerheblichen Vorteil gegenüber der linearen Regression (mit Encoding). Die Feature-Importances der beiden Kategorien sind vergleichsweise ähnlich, es kann also nicht direkt auf eine von beiden verzichtet werden. Insgesamt bleibt auch dieses Modell bzw. beide Modelle mit nur kategoriellen Features weit hinter den Modellen mit numerischen Features zurück.
del x_train
del x_val
del y_train
del y_val
Aufgabe A-7: Benchmark-Modelle mit numerischen und kategoriellen Features erstellen¶
Nach den bisher betrachteten Modellen, die entweder rein auf den numerischen bzw. kategoriellen Merkmalen bestanden, sollen in diesem Abschnitt kombinierte Modelle betrachtet werden.
Aufgabe A-7 a)¶
Analog zu den Aufgaben A-4 und A-6 ist ein lineares Modell mit Berücksichtigung sowohl der numerischen Merkmale als auch der kategoriellen Merkmale ocean_proximity und county_name zu erstellen. Die Ergebnisse der Modellierung sind auszugeben und zu kommentieren sowie die Ergebniswerte samt Laufzeit zu speichern. Wie ist dieses Modell im Vergleich zu den beiden anderen (linearen) zu bewerten?
Lineare Regression der Hauswerte mit sowohl numerischen als auch kategoriellen Merkmalen. Da lineare Regression ein Encoding erfordert, wird ein One-Hot-Encoding angewendet.
encoded_data=pd.get_dummies(data, columns=categorical_features, drop_first=True)
columns_removed=data.columns.difference(encoded_data.columns)
columns_added=encoded_data.columns.difference(data.columns)
data_train = encoded_data[encoded_data['sample'] == 'A']
x_train = data_train.drop(['median_house_value','sample'], axis=1)
y_train = data_train['median_house_value']
del data_train
data_val = encoded_data[encoded_data['sample'] == 'B']
x_val = data_val.drop(['median_house_value','sample'], axis=1)
y_val = data_val['median_house_value']
del data_val
del encoded_data
del columns_added, columns_removed
model_LR = sklearn.linear_model.LinearRegression()
starttime = time.time()
model_LR.fit(x_train, y_train)
fittime = time.time() - starttime
times['Complete Linear'] = fittime
save_times(times)
y_pred_LR = model_LR.predict(x_val)
print("Laufzeit:",fittime,"s")
del starttime
del fittime
Laufzeit: 0.11070752143859863 s
print("RMSE auf den Trainingsdaten (zur Erkennung von Overfitting):",rmse(y_train, model_LR.predict(x_train)))
RMSE auf den Trainingsdaten (zur Erkennung von Overfitting): 66285.48003351722
print("Ergebnisse auf den Validierungsdaten:")
print_results(y_val, y_pred_LR)
rmses['Complete Linear'] = rmse(y_val, y_pred_LR)
del y_pred_LR
del model_LR
Ergebnisse auf den Validierungsdaten: RMSE: 69098.04473488631 RMST (zum Vergleich): 114008.14636807823 R^2: 0.6326671677348177
print("Vergleich der RMSEs auf den Valdiierungsdaten")
display(rmses)
Vergleich der RMSEs auf den Valdiierungsdaten
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631}
Mit dem kompletten linearen Modell ist eine Verbesserung gegenüber den rein kategoriellen und numerischen linearen Modellen zu verzeichnen, es ist jedoch keine besonders große Verbesserung gegenüber dem numerischen Modell. Insbesondere ist das (naive) lineare Modell mit allen Merkmalen immer noch deutlich schlechter als das rein numerische Catboost. Da der RMSE auf den Trainingsdaten sogar niedriger ist als auf den Validierungsdaten scheint es sich hierbei nicht um Overfitting zu handeln.
del x_train
del x_val
del y_train
del y_val
Aufgabe A-7 b)¶
Analog zu den Aufgaben A-4 und A-6 ist CatBoost (ohne weitere Anpassungen) mit Berücksichtigung sowohl der numerischen Merkmale als auch der kategoriellen Merkmale ocean_proximity und county_name anzuwenden. Ein Ergebnisbericht ist auszugeben und die Ergebniswerte samt Laufzeit zu speichern. Zusätzlich ist eine grafische Feature Importance auszugeben und zusammen mit den Ergebniswerten zu kommentieren. Wie ist dieses Modell im Vergleich zu den beiden anderen (CatBoost) zu bewerten?
CatBoost-Modell mit allen Merkmalen. Hier wird absichtlich für CatBoost ein Datensatz ohne Encoding verwendet, weil das manuelle Encoding in Aufgabenteil B geschehen soll.
data_train = data[data['sample'] == 'A']
x_train = data_train.drop(['median_house_value','sample'], axis=1)
y_train = data_train['median_house_value']
del data_train
data_val = data[data['sample'] == 'B']
x_val = data_val.drop(['median_house_value','sample'], axis=1)
y_val = data_val['median_house_value']
del data_val
data_test = data[data['sample'] == 'C']
x_test = data_test.drop(['median_house_value','sample'], axis=1)
y_test = data_test['median_house_value']
del data_test
model_CB = catboost.CatBoostRegressor(verbose=False, cat_features=categorical_features, random_state=seed)
starttime = time.time()
model_CB.fit(x_train, y_train)
fittime = time.time() - starttime
times['Complete CatBoost'] = fittime
save_times(times)
y_pred_CB = model_CB.predict(x_val)
print("Laufzeit:",fittime,"s")
del starttime
del fittime
Laufzeit: 12.42766809463501 s
print("RMSE auf den Trainingsdaten (zur Erkennung von Overfitting):",rmse(y_train, model_CB.predict(x_train)))
RMSE auf den Trainingsdaten (zur Erkennung von Overfitting): 38038.35963547316
y_pred_CB = model_CB.predict(x_val)
print_results(y_val, y_pred_CB)
rmses['Complete Catboost'] = rmse(y_val, y_pred_CB)
del y_pred_CB
display(model_CB.get_params())
RMSE: 48361.638797917374 RMST (zum Vergleich): 114008.14636807823 R^2: 0.8200589506135714
{'loss_function': 'RMSE',
'verbose': False,
'random_state': 43,
'cat_features': ['ocean_proximity', 'county_name']}
print("Vergleich der RMSEs auf den Validierungsdaten")
display(rmses)
Vergleich der RMSEs auf den Validierungsdaten
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374}
Catboost mit allen Merkmalen ist deutlich besser als die lineare Regression, konnte allerdings die Performance gegenüber dem Catboost mit nur numerischen Merkmalen kaum verbessern. Ein Vergleich mit dem RMSE auf den Trainingsdaten, welcher etwas niedriger ist, lässt schließen, dass hier zumindest etwas Overfitting vorliegt.
# Auswertung auf den Testdaten für Aufgabe B-5 a)
y_pred_CB = model_CB.predict(x_test)
print_results(y_test, y_pred_CB)
rmsestest['Complete Catboost'] = rmse(y_test, y_pred_CB)
del y_pred_CB
RMSE: 45427.075729708304 RMST (zum Vergleich): 116233.22814781606 R^2: 0.8472543107282164
Feature Improtances:
importances = pd.DataFrame({'features': x_train.columns, 'importances': model_CB.get_feature_importance()}).sort_values(by='importances',ascending=False)
plt.figure(figsize=(16, 18))
sns.barplot(x='importances', y='features', data=importances)
plt.xlabel('Importance')
plt.ylabel('Merkmal')
plt.title('Feature Importance')
plt.show()
del importances
del model_CB
Anhand der Feature Importance kann man erkennen, warum das Modell nur eine geringe Verbesserung gegenüber dem CatBoost mit nur numerischen Merkmalen erzielt: Die wichtigsten Merkmale sind mit dem Medianeinkommen und den geografischen Standorten ohnehin numerisch, während die Ozeannähe auf dem vierten Platz immerhin noch eine ähnliche Signifikanz aufweisen kann.
del x_train
del x_val
del y_train
del y_val
Aufgabe A-7 c)¶
Wie sind die Laufzeiten und Prognosefehler der in dieser Aufgabe erstellten Modelle im Vergleich und wie sind die Werte und Unterschiede zu bewerten? Basierend auf den Ergebnissen dieses Abschnitts soll ein abschließendes Fazit gezogen werden.
Prognosefehler der Modelle via RMSE:
rmses
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374}
Laufzeitvergleich zwischen linearer Regression und CatBoost in Sekunden:
times
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501}
Als abschließendes Fazit kann man erkennen dass CatBoost gegenüber der linearne Regression in diesem Fall eine deutlich bessere Vorhersagegüte bietet, allerdings auch eine deutlich erhöhte Laufzeit erfordert.
Die Prognosegüte der linearen Regression ist vermutlich daher deutlich schlechter, weil sich vor allem der geografische Ort mittels der Längen- und Breitengrade nicht besonder gut in einem linearen Modell abbilden lässt (zumindest ohne komplexere Transformationen). CatBoost hat dieses Problem bei Entscheidungsbäumen naturgemäß nicht.
Der Unterschied in der Laufzeit ist (zumindest bei CPU-Durchführung) mehr als eine komplette Größenordnung und könnte bei einem entsprechend großen Datensatz zusammen mit schlechter Hardware oder Echtzeitanfordeurngen (beim Fitting) ein Ausschlusskriterium für CatBoost darstellen.
Teil B: Prognosemodelle optimieren und Overfitting verhindern¶
Aufgabe B-1: Encoding, Skalierung und Vergleich¶
Aufgabe B-1 a)¶
Die numerischen Daten sind geeignet zu skalieren. Das ausgewählte Skalierungsverfahren ist kurz zu beschreiben und die Auswahl kurz zu begründen. Es sind vor und nach Skalierung die gleichen sieben zufälligen Zeilen der Daten anzuzeigen (die gleichen Datensätze) und das Ergebnis ist zu bewerten.
Es wird MimMax-Scaling verwendet um die numerischen Daten zu skalieren.MinMax-Scaling skaliert die Differenz zwischen Minimum und Maximum einer numerischen Spalte auf 1 und verschiebt die Werte, sodass das Minimum bei 0 und das Maximum bei 1 liegt. Diese Methode wurde hauptsächlich gegenüber einer Standardisierung ausgewählt, weil es insbesondere für die weitere Anwendung in neuronalen Netzen effektiv und einfach zu handhaben ist.
sdata = data.copy()
featurescaler = sklearn.preprocessing.MinMaxScaler()
targetscaler = sklearn.preprocessing.MinMaxScaler()
numeric_features = numeric_columns.copy()
target_column = 'median_house_value'
numeric_features.remove(target_column)
sdata[numeric_features] = featurescaler.fit_transform(sdata[numeric_features])
sdata[target_column] = targetscaler.fit_transform(sdata[target_column].values.reshape(-1,1))
Ausgabe von 7 Zeilen vor und nach der Transformation
sample_index = sdata.sample(7, random_state=seed).index
display(data.loc[sample_index])
display(sdata.loc[sample_index])
| latitude | longitude | housing_median_age | population | median_income | median_house_value | ocean_proximity | county_name | sample | rooms_per_house | |
|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 34.05 | -117.27 | 34 | 849 | 3.1607 | 138200 | INLAND | San Bernardino | A | 4.743733 |
| 4402 | 33.88 | -117.57 | 35 | 1453 | 2.3160 | 119400 | INLAND | Riverside | C | 4.100467 |
| 1929 | 33.60 | -117.90 | 25 | 906 | 3.6538 | 500001 | <1H OCEAN | Orange | A | 5.222458 |
| 11551 | 35.37 | -118.97 | 34 | 1156 | 1.7197 | 48900 | INLAND | Kern | A | 4.377778 |
| 9882 | 34.20 | -118.36 | 14 | 1874 | 2.5267 | 231800 | <1H OCEAN | Los Angeles | C | 3.359571 |
| 3896 | 33.85 | -118.39 | 24 | 1830 | 5.7128 | 366200 | <1H OCEAN | Los Angeles | C | 5.307039 |
| 16638 | 37.83 | -122.25 | 52 | 534 | 3.1429 | 231600 | NEAR BAY | Alameda | B | 4.395189 |
| latitude | longitude | housing_median_age | population | median_income | median_house_value | ocean_proximity | county_name | sample | rooms_per_house | |
|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 0.160468 | 0.705179 | 0.647059 | 0.023711 | 0.183501 | 0.254022 | INLAND | San Bernardino | A | 0.027630 |
| 4402 | 0.142402 | 0.675299 | 0.666667 | 0.040640 | 0.125247 | 0.215259 | INLAND | Riverside | C | 0.023070 |
| 1929 | 0.112646 | 0.642430 | 0.470588 | 0.025309 | 0.217507 | 1.000000 | <1H OCEAN | Orange | A | 0.031024 |
| 11551 | 0.300744 | 0.535857 | 0.647059 | 0.032316 | 0.084123 | 0.069899 | INLAND | Kern | A | 0.025036 |
| 9882 | 0.176408 | 0.596614 | 0.254902 | 0.052440 | 0.139777 | 0.447011 | <1H OCEAN | Los Angeles | C | 0.017818 |
| 3896 | 0.139214 | 0.593625 | 0.450980 | 0.051207 | 0.359505 | 0.724123 | <1H OCEAN | Los Angeles | C | 0.031623 |
| 16638 | 0.562168 | 0.209163 | 1.000000 | 0.014883 | 0.182273 | 0.446598 | NEAR BAY | Alameda | B | 0.025159 |
Die besonders großen Werte einzelner Datensätze führen dazu, dass im Median die Werte von z.B. population und median_income nach der Skalierung in vielen Fällen sehr geringe Werte annimmt. Ansonsten ist die Transformation aber wie erwartet und (mit dem Wissen über obere Grenzwerte) ohne große Auffälligkeit.
Aufgabe B-1 b)¶
Die kategoriellen Daten sind geeignet in numerische Werte zu überführen. Das ausgewählte Encoding ist kurz zu beschreiben und die Auswahl kurz zu begründen. Für die sieben Zeilen aus Teilaufgabe a) sind die bearbeiteten Datensätze anzuzeigen und das Ergebnis ist zu bewerten. Es ist darauf einzugehen, wie viele Merkmale ggf. neu hinzugekommen sind und ob die Daten zu einem überbestimmten linearen Modell führen würden.
Wir verwenden One-Hot-Encoding um für jede Ausprägung jeweils ein boolsches Merkmal zu haben, welches anzeigt, ob das Merkmal zutrifft. Ausnahme bildet die erste Ausprägung, welche, um Redundanz zu vermeiden, kein eigenes Merkmal bekommt, sondern durch nicht-zutreffen aller anderen Merkmale codiert ist.
Die boolschen Merkmale werden in einer weiteren numerischen Verarbeitung typischerweise mit false=0 und true=1 verarbeitet.
Da es bei den Lankreisen keinen ordinalen Zusammenhang gibt ist One-Hot-Encoding einfach und effektiv um die Kategorien ohne weiteren Kontext zu repräsentieren.
esdata=pd.get_dummies(sdata, columns=categorical_features, drop_first=True)
columns_removed=sdata.columns.difference(esdata.columns)
columns_added=esdata.columns.difference(sdata.columns)
Darstellung der 7 Beispielzeilen vor und nach dem Encoding
display(sdata.loc[sample_index])
display(esdata.loc[sample_index])
del sample_index
| latitude | longitude | housing_median_age | population | median_income | median_house_value | ocean_proximity | county_name | sample | rooms_per_house | |
|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 0.160468 | 0.705179 | 0.647059 | 0.023711 | 0.183501 | 0.254022 | INLAND | San Bernardino | A | 0.027630 |
| 4402 | 0.142402 | 0.675299 | 0.666667 | 0.040640 | 0.125247 | 0.215259 | INLAND | Riverside | C | 0.023070 |
| 1929 | 0.112646 | 0.642430 | 0.470588 | 0.025309 | 0.217507 | 1.000000 | <1H OCEAN | Orange | A | 0.031024 |
| 11551 | 0.300744 | 0.535857 | 0.647059 | 0.032316 | 0.084123 | 0.069899 | INLAND | Kern | A | 0.025036 |
| 9882 | 0.176408 | 0.596614 | 0.254902 | 0.052440 | 0.139777 | 0.447011 | <1H OCEAN | Los Angeles | C | 0.017818 |
| 3896 | 0.139214 | 0.593625 | 0.450980 | 0.051207 | 0.359505 | 0.724123 | <1H OCEAN | Los Angeles | C | 0.031623 |
| 16638 | 0.562168 | 0.209163 | 1.000000 | 0.014883 | 0.182273 | 0.446598 | NEAR BAY | Alameda | B | 0.025159 |
| latitude | longitude | housing_median_age | population | median_income | median_house_value | sample | rooms_per_house | ocean_proximity_INLAND | ocean_proximity_ISLAND | ... | county_name_Sonoma | county_name_Stanislaus | county_name_Sutter | county_name_Tehama | county_name_Trinity | county_name_Tulare | county_name_Tuolumne | county_name_Ventura | county_name_Yolo | county_name_Yuba | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 7310 | 0.160468 | 0.705179 | 0.647059 | 0.023711 | 0.183501 | 0.254022 | A | 0.027630 | True | False | ... | False | False | False | False | False | False | False | False | False | False |
| 4402 | 0.142402 | 0.675299 | 0.666667 | 0.040640 | 0.125247 | 0.215259 | C | 0.023070 | True | False | ... | False | False | False | False | False | False | False | False | False | False |
| 1929 | 0.112646 | 0.642430 | 0.470588 | 0.025309 | 0.217507 | 1.000000 | A | 0.031024 | False | False | ... | False | False | False | False | False | False | False | False | False | False |
| 11551 | 0.300744 | 0.535857 | 0.647059 | 0.032316 | 0.084123 | 0.069899 | A | 0.025036 | True | False | ... | False | False | False | False | False | False | False | False | False | False |
| 9882 | 0.176408 | 0.596614 | 0.254902 | 0.052440 | 0.139777 | 0.447011 | C | 0.017818 | False | False | ... | False | False | False | False | False | False | False | False | False | False |
| 3896 | 0.139214 | 0.593625 | 0.450980 | 0.051207 | 0.359505 | 0.724123 | C | 0.031623 | False | False | ... | False | False | False | False | False | False | False | False | False | False |
| 16638 | 0.562168 | 0.209163 | 1.000000 | 0.014883 | 0.182273 | 0.446598 | B | 0.025159 | False | False | ... | False | False | False | False | False | False | False | False | False | False |
7 rows × 69 columns
print("Es wurden",columns_removed.size,"kategorielle Spalten entfernt.")
print("Durch das Encoding kamen",columns_added.size,"binäre Spalten hinzu.")
del columns_removed, columns_added
Es wurden 2 kategorielle Spalten entfernt. Durch das Encoding kamen 61 binäre Spalten hinzu.
Da wir beim Encoding die erste Spalte entfernt haben würde es hier kein überbestimmtes lineares Modell gebem. Aufgrund der vielen Ausprägungen der Landkreise gibt es viele neue Spalten durch das Encoding, die Dimension ist allerdings noch nicht in einem problematischen Bereich.
Aufgabe B-1 c)¶
Im Anschluss ist der nun skalierte und encodierte Datensatz in Trainingsdaten (sample = A) und Validierungsdaten (sample = B) aufzusplitten und die Zielvariable von den beschreibenden Variablen zu trennen.
Einordnung in Test- und Validierungsdatensatz
data_train = esdata[esdata['sample'] == 'A']
x_train = data_train.drop([target_column,'sample'], axis=1)
y_train = data_train[target_column]
del data_train
data_val = esdata[esdata['sample'] == 'B']
x_val = data_val.drop([target_column,'sample'], axis=1)
y_val = data_val[target_column]
del data_val
#für Aufgabe B-5
data_test = esdata[esdata['sample'] == 'C']
x_test = data_test.drop([target_column,'sample'], axis=1)
y_test = data_test[target_column]
del data_test
Aufgabe B-1 d)¶
An die so aufbereiteten Daten ist erneut ein CatBoost-Modell anzupassen. Dabei ist eine Zeitmessung durchzuführen und die Feature Importance auszugeben. Die Güte des Modells (RMSE) ist zunächst für den Trainingsdatensatz, im Anschluss für den Validierungsdatensatz zu ermitteln und die Ergebnisse sind zu vergleichen und zu interpretieren. Abschließend ist ein Vergleich mit dem CatBoost-Modell aus Aufgabe A-7 vorzunehmen und auf Laufzeit, Modellgüte (Validierung) und ggf. Auffälligkeiten im Vergleich der Feature Importances einzugehen, mögliche Ursachen zu ergründen und eine Bewertung vorzunehmen.
CatBoost an die skalierten und encodierten Daten anpassen
model_CB = catboost.CatBoostRegressor(random_state=seed,verbose=False)
starttime = time.time()
model_CB.fit(x_train, y_train)
fittime = time.time() - starttime
times['Scaled Encoded CatBoost'] = fittime
save_times(times)
y_val_pred_CB = model_CB.predict(x_val)
print("Laufzeit:",fittime,"s")
del starttime
del fittime
Laufzeit: 6.519219160079956 s
Um die RMSEs mit den vorherigen Modellen vergleichbar zu machen, müssen die Zielwerte wieder zurückskaliert werden.
y_train_reformed = targetscaler.inverse_transform(y_train.values.reshape(-1,1))
y_train_pred_CB = model_CB.predict(x_train)
y_train_pred_reformed = targetscaler.inverse_transform(y_train_pred_CB.reshape(-1,1))
print("RMSE auf den Trainingsdaten (zur Erkennung von Overfitting):",rmse(y_train_reformed, y_train_pred_reformed))
RMSE auf den Trainingsdaten (zur Erkennung von Overfitting): 37573.54660791557
y_val_pred_CB_reformed = targetscaler.inverse_transform(y_val_pred_CB.reshape(-1,1))
y_val_reformed = targetscaler.inverse_transform(y_val.values.reshape(-1,1))
print_results(y_val_reformed, y_val_pred_CB_reformed)
rmses['Scaled Encoded Catboost'] = rmse(y_val_reformed, y_val_pred_CB_reformed)
del y_val_pred_CB_reformed
del y_val_reformed
display(model_CB.get_params())
RMSE: 48406.630856312855 RMST (zum Vergleich): 114008.14636807823 R^2: 0.8197239874352333
{'loss_function': 'RMSE', 'verbose': False, 'random_state': 43}
#für Aufgabe B-5
y_test_pred_CB = model_CB.predict(x_test)
y_test_pred_CB_reformed = targetscaler.inverse_transform(y_test_pred_CB.reshape(-1,1))
y_test_reformed = targetscaler.inverse_transform(y_test.values.reshape(-1,1))
print_results(y_test_reformed, y_test_pred_CB_reformed)
rmsestest['Scaled Encoded Catboost'] = rmse(y_test_reformed, y_test_pred_CB_reformed)
del y_test_pred_CB, y_test_pred_CB_reformed
del y_test_reformed
RMSE: 46154.39938978077 RMST (zum Vergleich): 116233.22814781606 R^2: 0.8423239947821303
importances = pd.DataFrame({'features': x_train.columns, 'importances': model_CB.get_feature_importance()}).sort_values(by='importances',ascending=False)
plt.figure(figsize=(16, 18))
sns.barplot(x='importances', y='features', data=importances)
plt.xlabel('Importance')
plt.ylabel('Merkmal')
plt.title('Feature Importance')
plt.show()
del importances
del model_CB
print("Fitting-Laufzeiten")
times
Fitting-Laufzeiten
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956}
print("Vergleich der RMSEs auf den Validierungsdaten")
display(rmses)
Vergleich der RMSEs auf den Validierungsdaten
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855}
Insgesamt wir durch die Skalierung und das Encoding im eigentlichen Fitting Laufzeit gespart, die tatsächliche Prognosegüte ist aber nicht besser geworden. Grund dafür ist, dass Catboost von sich aus bereits mit kategoriellen Daten umgehen kann -- da unser One-Hot-Encoding nicht so optimiert wie das interne Encoding von CatBoost ist, erhöhen wir damit die Dimensionalität ohne Mehrwert. Hier ist, auch für spätere Interpretationen der numerischen Werte, das CatBoost auf unskalierten und unkodierten Daten zu bevorzugen.
Auffällig ist, dass die Feature Importance der Ozeannähe aus Aufgabe A-7 scheinbar fast vollständig durch die Ausprägung INLAND und in etwas geringerem Maße NEAR_OCEAN dominiert wird. Nach weiteren numerischen Merkmalen kommen anschließend mit wenigen Landkreisen dann wieder kategorische Merkmale in der Rangfolge. Insgesamt bleibt das Modell durch numerische Merkmale dominiert und gibt nur gewissen Ausprägungen der kategorischen Merkmale eine relative Wichtigkeit, da nur diese eine starke Aussagekraft zu besitzen scheinen.
del y_train_pred_CB
del y_train_pred_reformed
del y_train_reformed
del x_train
del x_val
del y_train
del y_val
del y_val_pred_CB
del targetscaler, featurescaler
Aufgabe B-2: Under-/Overfitting am Beispiel von Neuronalen Netzen¶
Aufgabe B-2 a)¶
Es ist ein Teildatensatz aus dem oben erstellten encodierten und skalierten Datensatz zu erstellen, der nur die ersten 100 Zeilen des Datensatzes enthält. Damit sind 3 einfache Neuronale Netze unterschiedlicher Größe mit Keras/Tensorflow sequenziell zu definieren, zu kompilieren, anzuzeigen (via summary) und jeweils über eine hohe Anzahl von Epochen (> 2000) zu trainieren: • klein: 1.000 bis 10.000 Gewichte • mittel: 50.000 bis 100.000 Gewichte • groß: mindestens 1.000.000 Gewichte Die history der einzelnen Modelle ist in einem (einzigen) Diagramm zu plotten. Dabei soll der Verlauf des RMSE auf den Trainingsdaten sowie der Verlauf des RMSE auf der internen Modellvalidierung (Split 0.2) dargestellt werden. Die Ergebnisse sind zu beschreiben und zu bewerten.
#Zurücksetzen der Seeds
random.seed(seed)
np.random.seed(seed)
tf.random.set_seed(seed)
Erstellung des Datensatzes. Da die vorhandenen Daten nach Lokalität geordnet sind und der interne Training-Validierungs-Split in Keras keine weitere Mischung vornimmt (sondern ab einem gewissen Punkt abschneidet), müssen die Trainingsdaten hier in eine zufällige Reihenfolge gebracht werden, um dem Bias in der Location entgegenzuwirken.
Es ist zu erwähnen, dass die Forderung, die ersten 100 Datensätze zu nehmen, bereits einen sehr großen geografischen Bias beinhaltet und die so erzeugten Modelle nicht gut auf den restlichen Datensatz verallgemeinern.
first_esdata_randomized=esdata[:100].sample(frac=1, random_state=seed)
y_train = first_esdata_randomized['median_house_value']
x_train = first_esdata_randomized.drop(['sample','median_house_value'], axis=1)
del first_esdata_randomized
Wir nutzen eine einfache Heuristik für die Netzwerkstruktur: DIe Anzahl der Layer skaliert logarithmisch mit der Anzahl der Neuronen. Die Anzahl der Neuronen pro Layer ist konstant, mit einem Ergebnislayer welches ein einziges Neuron umfasst.
Als Verlustmetrik nehmen wir direkt den MSE um daraus einfach den RMSE zurückgewinnen zu können.
Aufbau und Struktur der drei Netzwerke:
print("Erstes Netzwerk (klein)")
model1=simple_network(100, num_features=x_train.columns.size)
model1.summary()
print("")
print("Zweites Netzwerk (mittel)")
model2=simple_network(650, num_features=x_train.columns.size)
model2.summary()
print("")
print("Drittes Netzwerk (groß)")
model3=simple_network(3000, num_features=x_train.columns.size)
model3.summary()
Erstes Netzwerk (klein)
Model: "sequential"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
dense (Dense) (None, 25) 1700
dense_1 (Dense) (None, 25) 650
dense_2 (Dense) (None, 25) 650
dense_3 (Dense) (None, 25) 650
dense_4 (Dense) (None, 1) 26
=================================================================
Total params: 3676 (14.36 KB)
Trainable params: 3676 (14.36 KB)
Non-trainable params: 0 (0.00 Byte)
_________________________________________________________________
Zweites Netzwerk (mittel)
Model: "sequential_1"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
dense_5 (Dense) (None, 108) 7344
dense_6 (Dense) (None, 108) 11772
dense_7 (Dense) (None, 108) 11772
dense_8 (Dense) (None, 108) 11772
dense_9 (Dense) (None, 108) 11772
dense_10 (Dense) (None, 108) 11772
dense_11 (Dense) (None, 1) 109
=================================================================
Total params: 66313 (259.04 KB)
Trainable params: 66313 (259.04 KB)
Non-trainable params: 0 (0.00 Byte)
_________________________________________________________________
Drittes Netzwerk (groß)
Model: "sequential_2"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
dense_12 (Dense) (None, 375) 25500
dense_13 (Dense) (None, 375) 141000
dense_14 (Dense) (None, 375) 141000
dense_15 (Dense) (None, 375) 141000
dense_16 (Dense) (None, 375) 141000
dense_17 (Dense) (None, 375) 141000
dense_18 (Dense) (None, 375) 141000
dense_19 (Dense) (None, 375) 141000
dense_20 (Dense) (None, 1) 376
=================================================================
Total params: 1012876 (3.86 MB)
Trainable params: 1012876 (3.86 MB)
Non-trainable params: 0 (0.00 Byte)
_________________________________________________________________
Zum Trainieren werden hier 2024 Epochen genutzt.
epochs = 2024
trained_model_1, history_1 = train_model(model1, x_train, y_train, epochs, 'NN_klein',retrain_models, True, nn_verbosity)
/usr/local/lib/python3.10/dist-packages/keras/src/engine/training.py:3103: UserWarning: You are saving your model as an HDF5 file via `model.save()`. This file format is considered legacy. We recommend using instead the native Keras format, e.g. `model.save('my_model.keras')`.
saving_api.save_model(
trained_model_2, history_2 = train_model(model2, x_train, y_train, epochs,'NN_mittel',retrain_models, True, nn_verbosity)
trained_model_3, history_3 = train_model(model3, x_train, y_train, epochs,'NN_groß',retrain_models, True, nn_verbosity)
Um die orginale Skalierung des RMSE vor der MinMax-Skalierung zurückzugewinnen benötigen wir den Faktor $(y_{max} - y_{min})^2$ der unskalierten Werte, denn $$MSE_{org}=(y_{max} - y_{min})^2 MSE_{skaliert}$$ $$RMSE_{org}=\sqrt{MSE_{org}}$$
MSE_scaling_factor = (data['median_house_value'].max() - data['median_house_value'].min()) **2
scaled_loss_1 = [math.sqrt(val * MSE_scaling_factor) for val in history_1.history['mean_squared_error'] ]
scaled_val_loss_1 = [math.sqrt(val * MSE_scaling_factor) for val in history_1.history['val_mean_squared_error'] ]
scaled_loss_2 = [math.sqrt(val * MSE_scaling_factor) for val in history_2.history['mean_squared_error' ] ]
scaled_val_loss_2 = [math.sqrt(val * MSE_scaling_factor) for val in history_2.history['val_mean_squared_error'] ]
scaled_loss_3 = [math.sqrt(val * MSE_scaling_factor) for val in history_3.history['mean_squared_error'] ]
scaled_val_loss_3 = [math.sqrt(val * MSE_scaling_factor) for val in history_3.history['val_mean_squared_error'] ]
plt.figure(figsize=(20, 10))
plt.plot(scaled_loss_1,linestyle='--',color='lightblue')
plt.plot(scaled_val_loss_1,linestyle='-',color='darkblue')
plt.plot(scaled_loss_2,linestyle='--',color='lightgreen')
plt.plot(scaled_val_loss_2,linestyle='-',color='darkgreen')
plt.plot(scaled_loss_3,linestyle='--',color='red')
plt.plot(scaled_val_loss_3,linestyle='-',color='darkred')
plt.title('RMSE nach Epoche')
plt.xlabel('Epoche')
plt.ylabel('RMSE')
plt.ylim(0,50000)
plt.legend(['Training Netzwerk 1 (klein)', 'Validierung Netzwerk 1 (klein)',
'Training Netzwerk 2 (mittel)', 'Validierung Netzwerk 2 (mittel)',
'Training Netzwerk 3 (groß)', 'Validierung Netzwerk 3 (groß)'],
loc='upper left')
plt.show()
Hier sind die drei Netzwerke (klein=blau, mittel=grün und groß=rot) jeweils mit ihrem zurückskalierten Trainings-RMSE (hell) und dem Validierungs-RMSE (dunkel) dargestellt. Es wird keinerlei Regularisierung verwendet.
Es ist auffällig, dass das kleine Netzwrk auch nach einer großen Anzahl an Epochen bei einem stabilen RMSE verbleibt. Bei den mittleren und großen Netzwerken ist schnell ein Effekt von Overfitting zu sehen, bei dem der RMSE auf den Trainingsdaten verbessert wird, auf den Validierungsdaten dafür aber deutlich schlechter wird.
Das kleine Netzwerk performt auch deswegen besser, weil durch den sehr kleinen Datensatz gar nicht mehr generalisierbare Informationen in den Trainingsdaten vorhanden sind, insbesondere im Vergleich mit den zu trainierenden Parametern.
Da wir nur die 100 ersten Datensätze verwendet haben und die Datensätze geografisch sortiert sind haben wir hier auch mehr uniformität als im gesamten Datensatz, womit auch ein geringerer RMSE einher geht. Bei Verwendung des vollen Datensatzes ist zu erwarten, dass der RMSE entsprechend höher ausfallen wird.
print("Anzahl Informationswerte im Trainingssatz: ",0.8 * x_train.columns.size * len(x_train))
Anzahl Informationswerte im Trainingssatz: 5360.0
del model1
del model2
del model3
del trained_model_1
del trained_model_2
del trained_model_3
del scaled_loss_1
del scaled_val_loss_1
del scaled_loss_2
del scaled_val_loss_2
del scaled_loss_3
del scaled_val_loss_3
del history_1
del history_2
del history_3
del x_train
del y_train
del MSE_scaling_factor
Aufgabe B-2 b)¶
Es soll ein Neuronales Netz für den kompletten Trainingsdatensatz (sample = A) erstellt werden. Hierbei sind bei der Architektur des Netzes ggf. Erkenntnisse aus Teilaufgabe a) miteinzubeziehen und die vorgenommenen Änderungen an der Architektur des Netzes zu beschreiben. Die Größe und Architektur des Neuronalen Netzes muss ausreichend groß sein, um die Komplexität des Datensatzes abbilden zu können. Overfitting soll durch geeignete Regularisierung verhindert werden. Die Daten sind zu fitten und die Lernkurve aus (interner) Trainings- und Modellvalidierung ist zu plotten. Im Anschluss ist das Modell mit den „externen“ Validierungsdaten (sample = B) zu validieren. Falls die Prognosegüte dieses hochparametrigen Modells schlechter als die des einfachen linearen Modells aus Aufgabe A-7 a) ist, sind die Ursachen zu suchen und zu beheben. Dieses Ergebnis (RMSE) ist samt Laufzeit anzuzeigen, dem Gesamtvergleich hinzuzufügen und zu bewerten.
#Zurücksetzen der Seeds
random.seed(seed)
np.random.seed(seed)
tf.random.set_seed(seed)
Aufbau der Trainingsdaten. Wir oben randomisieren wir die Reihenfolge der Trainingsdaten für eine bessere (interne) Validierung.
resdata=esdata.sample(frac=1, random_state=seed)
x_train = resdata[resdata['sample']=='A'].drop(['sample', target_column], axis=1)
y_train=resdata[resdata['sample']=='A'][target_column]
Eine grobe Einschätzung, wie viel tatsächliche Informationen in unserem Trainingsdatensatz vorliegt (Anzahl Datenpunkte mal Anzahl Merkmale)
len(esdata[esdata['sample']=='A'])*(esdata.columns.size-1)
980628
Das Netzwerk wird dieses mal mit Parametern für die Learning-Rate und die Regularisierung definiert, ansonsten wie oben mit logarithmischer Anzahl der Layer und sinkenden Neuronen pro Layer.
Wir wählen eine relativ geringe Regularisierungsstärke, da diese einen großen Effekt auf den Loss besitzt. Außerdem kann man damit die weitere Strategie dargestellt werden, die wir gegen Overfitting verwenden: Early Stopping.
model = simple_regularized_network(middle_neurons=338, num_features=x_train.columns.size, reg_strength=0.0001, learning_rate=0.001)
model.summary()
Model: "sequential_3"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
dense_21 (Dense) (None, 67) 4556
dense_22 (Dense) (None, 67) 4556
dense_23 (Dense) (None, 67) 4556
dense_24 (Dense) (None, 67) 4556
dense_25 (Dense) (None, 67) 4556
dense_26 (Dense) (None, 1) 68
=================================================================
Total params: 22848 (89.25 KB)
Trainable params: 22848 (89.25 KB)
Non-trainable params: 0 (0.00 Byte)
_________________________________________________________________
Wir nutzen wieder >2000 Epochen, diese werden jedoch nicht erreicht. Durch Early-Stopping wird das Training abgebrochen, sobald sich der Loss auf den internen Validierungsdaten über eine Epochendauer nicht verbessert, wonach anschließend die besten gefundenen Gewichte verwendet werden.
Das Modell wird trainiert und die Laufzeit gespeichert.
starttime = time.time()
trained_model, history = train_model_early_stopping(model=model, x_train=x_train, y_train=y_train, epochs=2024, name='NN_data',retrain_model=retrain_models, verbosity=nn_verbosity)
fittime = time.time() - starttime
if retrain_models:
times['Complete Neural Network'] = fittime
save_times(times)
if 'Complete Neural Network' in times:
print("Laufzeit:",times['Complete Neural Network'],"s")
/usr/local/lib/python3.10/dist-packages/keras/src/engine/training.py:3103: UserWarning: You are saving your model as an HDF5 file via `model.save()`. This file format is considered legacy. We recommend using instead the native Keras format, e.g. `model.save('my_model.keras')`.
saving_api.save_model(
Laufzeit: 150.87731003761292 s
Plotten der History
MSE_scaling_factor = (data['median_house_value'].max() - data['median_house_value'].min()) **2
scaled_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['mean_squared_error'] ]
scaled_val_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['val_mean_squared_error'] ]
plt.figure(figsize=(20, 10))
plt.plot(scaled_loss,linestyle='--',color='lightblue')
plt.plot(scaled_val_loss,linestyle='-',color='darkblue')
plt.title('RMSE nach Epoche')
plt.xlabel('Epoche')
plt.ylabel('RMSE')
plt.grid(True)
plt.ylim(40000,70000)
plt.legend(['Training', 'Validierung'],
loc='upper left')
plt.show()
display(min(scaled_val_loss))
54823.60573524983
display(min(scaled_loss))
49966.43565167835
Man sieht dass hier eindeutig der Fehler der Trainingsdaten sinkt, während sich der Fehler der internen Validiserungsdaten nach einiger Zeit nur noch unwesentlich verändert. Zu Beginn werden aber echte Effekte gefunden, da sich auch die Validierungsergebnisse leicht verbessern. Ab einem gewissen Punkt kann der Validierungsfehler nicht mehr verbessert werden, wodurch Early Stopping getriggert und das Training beendet wird.
Nun wird das Modell auf den externen Validierungsdaten getestet
x_val = resdata[resdata['sample']=='B'].drop(['sample', target_column], axis=1)
y_val=resdata[resdata['sample']=='B'][target_column]
MSE = trained_model.evaluate(x_val.values.astype(float), y_val.values.astype(float))[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmses['Complete Neural Network'] = RMSE
display(rmses['Complete Neural Network'])
del MSE
del RMSE
del x_train
del x_val
del y_train
del y_val
98/98 [==============================] - 0s 2ms/step - loss: 0.0160 - mean_squared_error: 0.0137
56852.606799089896
#für Aufgabe B-5
x_test = resdata[resdata['sample']=='C'].drop(['sample', target_column], axis=1)
y_test=resdata[resdata['sample']=='C'][target_column]
MSE = trained_model.evaluate(x_test.values.astype(float), y_test.values.astype(float))[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmsestest['Complete Neural Network'] = RMSE
display(rmsestest['Complete Neural Network'])
del MSE
del RMSE
del x_test
del y_test
98/98 [==============================] - 0s 2ms/step - loss: 0.0152 - mean_squared_error: 0.0129
55160.98627873368
Die anderen Modell-RMSEs und Laufzeiten auf den Validierungsdaten zum Vergleich
print("RMSEs:")
display(rmses)
print("Laufzeiten:")
display(times)
RMSEs:
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855,
'Complete Neural Network': 56852.606799089896}
Laufzeiten:
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956,
'Complete Neural Network': 150.87731003761292}
Die Modellperformance ist zwar besser als das einfache lineare Modell, bleibt aber deutlich hinter CatBoost zurück. Angesichts des erforderlichen Tunings und der hohen Laufzeitanforderung des Trainings für neuronale Netzwerke liefert das Neuronale Netz zwar einen validen Ansatz, ist in diesem Anwendungsfall alleridngs nicht das erste Mittel der Wahl.
del MSE_scaling_factor
del fittime
del model
del resdata
del scaled_loss
del scaled_val_loss
del starttime
del trained_model
Aufgabe B3: Neuronales Netz mit Embeddings¶
Aufgabe B-3 a)¶
Die kategoriellen Merkmale ocean_proximity und county_name sind über zwei-dimensionale Embeddings in ein Neuronales Netz einzufügen. Die verborgenen Schichten des Netzes sollen dabei aus dem Neuronalen Netz aus Aufgabe B-2 b) übernommen werden. Die nötigen Anpassungen der Modellumsetzung sind detailliert vorab zu beschreiben und dann durchzuführen. Die Trainingsdaten sind zu fitten und die Konvergenz ist zu bewerten. Falls Letztere nicht gegeben ist, so sind die notwendigen Veränderungen durchzuführen.
Wegen des Embeddings nutzen wir hier nur den skalierten Datensatz ohne One-Hot Encoding.
training_data=sdata[sdata['sample']=='A'].sample(frac=1.0, random_state=seed)
x_train=training_data.drop([target_column,'sample'],axis=1)
y_train=training_data[target_column]
del training_data
Wir verwenden die gleichen einfachen Netzwerkaufbau mit einer konstanten Anzahl an Neuronen pro innerer Schicht und einer Anzahl von Schichten die logarithmisch mit der Gesamtanzahl der Neuronen steigt. Da wir ebenfalls die gleiche Anzahl an inneren Neuronen verwenden ist der Aufbau somit identisch zur Aufgabe B-2 b).
middle_neurons=338
learning_rate=0.001
reg_strength=0.0001
layer_shape = simple_layer_shape(middle_neurons)
display(layer_shape)
[67, 67, 67, 67, 67]
Die bisherigen Inputs aus dem One-Hot-Encoding fallen weg. Statdessen werden zwei Input-Layer, jeweils für die Kategorien der Ozeannnähe und der Landkreise angelegt. Hinter diesen liegt jeweils ein 2-dimensionales Encoding-Layer. Jedes dieser Encoding-Layer muss wieder in ein Paar aus 1-dimensionalen Werten überführt werden, weshalb diese jeweils mit Flatten-Layer verbunden sind.
Nun liefern beide Flatten-Layer jeweils zwei Outputs und stellen die embeddeten Informationen der kategoriellen Variablen dar. Daneben wird das Input-Layer für die numerischen Features angelegt und mit einem Concatenation-Layer zu einem Gesamtlayer als Input für das dichte neuronale Netz zusammengefügt.
Das dichte neuronale Netz hat, wie oben beschrieben, den identischen einfachen Aufbau wie im vorherigen Aufgabenteil.
embedding_size = 2
cat_input_1 = keras.layers.Input(shape=(1,))
cat_input_2 = keras.layers.Input(shape=(1,))
embed_layer1 = keras.layers.Embedding(input_dim=x_train['ocean_proximity'].nunique(), output_dim=embedding_size)(cat_input_1)
embed_layer2 = keras.layers.Embedding(input_dim=x_train['county_name'].nunique(), output_dim=embedding_size)(cat_input_2)
flatten_layer_1 = keras.layers.Flatten()(embed_layer1)
flatten_layer_2 = keras.layers.Flatten()(embed_layer2)
num_input = keras.layers.Input(shape=(len(numeric_features),))
concatenated_inputs = keras.layers.Concatenate()([flatten_layer_1, flatten_layer_2, num_input])
previous_layer=concatenated_inputs
for nodes in layer_shape:
next_layer = keras.layers.Dense(nodes, activation='relu',kernel_regularizer=keras.regularizers.l2(reg_strength))(previous_layer)
previous_layer = next_layer
output_layer = keras.layers.Dense(1, activation='relu',kernel_regularizer=keras.regularizers.l2(reg_strength))(previous_layer)
model = keras.Model(inputs=[cat_input_1, cat_input_2, num_input], outputs=output_layer)
optimizer = keras.optimizers.Adam(learning_rate=learning_rate)
model.compile(optimizer=optimizer,
loss='mean_squared_error',
metrics=['mean_squared_error'])
del cat_input_1, cat_input_2
del embed_layer1, embed_layer2
del flatten_layer_1, flatten_layer_2
del num_input, concatenated_inputs
del previous_layer, next_layer, output_layer, optimizer
model.summary()
Model: "model"
__________________________________________________________________________________________________
Layer (type) Output Shape Param # Connected to
==================================================================================================
input_1 (InputLayer) [(None, 1)] 0 []
input_2 (InputLayer) [(None, 1)] 0 []
embedding (Embedding) (None, 1, 2) 10 ['input_1[0][0]']
embedding_1 (Embedding) (None, 1, 2) 116 ['input_2[0][0]']
flatten (Flatten) (None, 2) 0 ['embedding[0][0]']
flatten_1 (Flatten) (None, 2) 0 ['embedding_1[0][0]']
input_3 (InputLayer) [(None, 6)] 0 []
concatenate (Concatenate) (None, 10) 0 ['flatten[0][0]',
'flatten_1[0][0]',
'input_3[0][0]']
dense_27 (Dense) (None, 67) 737 ['concatenate[0][0]']
dense_28 (Dense) (None, 67) 4556 ['dense_27[0][0]']
dense_29 (Dense) (None, 67) 4556 ['dense_28[0][0]']
dense_30 (Dense) (None, 67) 4556 ['dense_29[0][0]']
dense_31 (Dense) (None, 67) 4556 ['dense_30[0][0]']
dense_32 (Dense) (None, 1) 68 ['dense_31[0][0]']
==================================================================================================
Total params: 19155 (74.82 KB)
Trainable params: 19155 (74.82 KB)
Non-trainable params: 0 (0.00 Byte)
__________________________________________________________________________________________________
starttime = time.time()
trained_model, history = train_embed_model(model,x_train,y_train,epochs=2024,name='NN_embed',retrain_model=retrain_models,verbosity=nn_verbosity)
fittime = time.time() - starttime
if retrain_models:
times['Embedding Neural Network'] = fittime
save_times(times)
if 'Embedding Neural Network' in times:
print("Laufzeit:",times['Embedding Neural Network'],"s")
del starttime
del fittime
/usr/local/lib/python3.10/dist-packages/keras/src/engine/training.py:3103: UserWarning: You are saving your model as an HDF5 file via `model.save()`. This file format is considered legacy. We recommend using instead the native Keras format, e.g. `model.save('my_model.keras')`.
saving_api.save_model(
Laufzeit: 203.7618260383606 s
MSE_scaling_factor = (data['median_house_value'].max() - data['median_house_value'].min()) **2
scaled_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['mean_squared_error'] ]
scaled_val_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['val_mean_squared_error'] ]
plt.figure(figsize=(20, 10))
plt.plot(scaled_loss,linestyle='--',color='lightblue')
plt.plot(scaled_val_loss,linestyle='-',color='darkblue')
plt.title('RMSE nach Epoche')
plt.xlabel('Epoche')
plt.ylabel('RMSE')
plt.grid(True)
plt.ylim(40000,70000)
plt.legend(['Training', 'Validierung'],
loc='upper left')
<matplotlib.legend.Legend at 0x7b40fca14a00>
Hier greift ebenfalls wieder das Early-Stopping, sobald der Validierungsfehler über viele Epochen hinweg nicht verbesser werden kann.
display(min(scaled_val_loss))
54335.74381974795
display(min(scaled_loss))
50749.03177217832
del embedding_size
del layer_shape
del middle_neurons
del model
del nodes
del reg_strength
del scaled_loss
del scaled_val_loss
Aufgabe B-3 b)¶
Der RMSE der Validierungsdaten des Netzes mit Embeddings ist samt Fitting-Laufzeit zu ermitteln und dem Modellvergleich hinzuzufügen. Das Neuronale Netz mit Embeddings ist mit dem Neuronalen Netz aus Aufgabe B-2 b) zu vergleichen und zu bewerten.
validation_data=sdata[data['sample']=='B'].sample(frac=1.0, random_state=seed)
x_val=validation_data.drop([target_column,'sample'],axis=1)
y_val=validation_data[target_column]
del validation_data
MSE = trained_model.evaluate([x_val['ocean_proximity'].cat.codes,x_val['county_name'].cat.codes, x_val[numeric_features].values.astype(float)], y_val.values.astype(float),verbose=0)[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmses['Embedding Neural Network'] = RMSE
print("RMSE auf Validierungsdaten: ",rmses['Embedding Neural Network'])
RMSE auf Validierungsdaten: 56651.739797472845
rmses
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855,
'Complete Neural Network': 56852.606799089896,
'Embedding Neural Network': 56651.739797472845}
times
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956,
'Complete Neural Network': 150.87731003761292,
'Embedding Neural Network': 203.7618260383606}
Das Embedding Neural Network kann keine signifikante Performance-Verbesserung gegenüber dem Neural Network mit One-Hot-Encoding erzielen. Es ist zu erwähnen, dass durch das nicht benötigte One-Hot-Encoding auch weniger Parameter bei den Inputs benötigt werden.
Bei der Laufzeit ist das das Neural Network mit Embedding etwas langsamer als das Neural Network mit One-Hot-Encoding, da das Encoding in jedem Schritt zusätzlich mittrainiert wird.
Falls der Encoding-Aufwand und die höhere Dimensionalität (inklusive Trainingsaufwand) kein Problem darstellen, lässt sich mit Encoding (ohne Embedding) eine gute Modellperformance erzielen. Beim Embedding hat man allerdings auch eine bessere Interpretierbarkeit durch die Embeddings, die ihrerseits in einem Anwendungsfall wichtig sein können.
del MSE
del RMSE
#für Aufgabe B-5
test_data=sdata[data['sample']=='C'].sample(frac=1.0, random_state=seed)
x_test=test_data.drop([target_column,'sample'],axis=1)
y_test=test_data[target_column]
MSE = trained_model.evaluate([x_test['ocean_proximity'].cat.codes,x_test['county_name'].cat.codes, x_test[numeric_features].values.astype(float)], y_test.values.astype(float),verbose=0)[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmsestest['Embedding Neural Network'] = RMSE
del MSE
del RMSE
del MSE_scaling_factor
del x_test, y_test
del test_data
Aufgabe B-3 c)¶
Die Embeddings der kategoriellen Merkmale sind aus dem Modell zu extrahieren und grafisch darzustellen.
feature_name='ocean_proximity'
ocean_prox_embedding = trained_model.layers[2]
ocean_prox_input = trained_model.input[0]
embedding_model_ocean = keras.models.Model(inputs=ocean_prox_input, outputs=ocean_prox_embedding.output)
plot_embedding(feature_name='ocean_proximity', embedding_model=embedding_model_ocean, x_train=x_train, title='Embeddings für die Ozeannähe')
1/1 [==============================] - 0s 149ms/step
feature_name='county_name'
county_name_embedding = trained_model.layers[3]
county_name_input = trained_model.input[1]
embedding_model_county = keras.models.Model(inputs=county_name_input, outputs=county_name_embedding.output)
plot_embedding(feature_name=feature_name, embedding_model=embedding_model_county, x_train=x_train, title='Embeddings für die Landkreise',x_size=30, y_size=20, fontsize=12)
2/2 [==============================] - 0s 6ms/step
del ocean_prox_embedding
del ocean_prox_input
del county_name_embedding
del county_name_input
Aufgabe B-3 d)¶
Das in Teilaufgabe a) aufgebaute Netz mit Embeddings soll als „clone“ nochmals, nun aber mit einem anderen Modellnamen und anderen (random) Startgewichten erzeugt, neu gefitted und die Konvergenz überprüft und ggf. hergestellt werden. Die entsprechenden Embedding-Gewichte sind zu extrahieren und in einer verbundenen Grafik den Embeddings aus Teilaufgabe c) gegenüberzustellen (z. B. 2x2 Grafiken über-/nebeneinander) und zu interpretieren.
model_clone = keras.models.clone_model(trained_model)
del trained_model
optimizer = keras.optimizers.Adam(learning_rate=learning_rate)
model_clone.compile(optimizer=optimizer,
loss='mean_squared_error',
metrics=['mean_squared_error'])
model_clone.summary()
Model: "model"
__________________________________________________________________________________________________
Layer (type) Output Shape Param # Connected to
==================================================================================================
input_1 (InputLayer) [(None, 1)] 0 []
input_2 (InputLayer) [(None, 1)] 0 []
embedding (Embedding) (None, 1, 2) 10 ['input_1[0][0]']
embedding_1 (Embedding) (None, 1, 2) 116 ['input_2[0][0]']
flatten (Flatten) (None, 2) 0 ['embedding[0][0]']
flatten_1 (Flatten) (None, 2) 0 ['embedding_1[0][0]']
input_3 (InputLayer) [(None, 6)] 0 []
concatenate (Concatenate) (None, 10) 0 ['flatten[0][0]',
'flatten_1[0][0]',
'input_3[0][0]']
dense_27 (Dense) (None, 67) 737 ['concatenate[0][0]']
dense_28 (Dense) (None, 67) 4556 ['dense_27[0][0]']
dense_29 (Dense) (None, 67) 4556 ['dense_28[0][0]']
dense_30 (Dense) (None, 67) 4556 ['dense_29[0][0]']
dense_31 (Dense) (None, 67) 4556 ['dense_30[0][0]']
dense_32 (Dense) (None, 1) 68 ['dense_31[0][0]']
==================================================================================================
Total params: 19155 (74.82 KB)
Trainable params: 19155 (74.82 KB)
Non-trainable params: 0 (0.00 Byte)
__________________________________________________________________________________________________
starttime = time.time()
cloned_trained_model, history = train_embed_model(model_clone,x_train,y_train,epochs=2024, name='NN_cloned',retrain_model=retrain_models,do_save_model=True, verbosity=nn_verbosity)
fittime = time.time() - starttime
if retrain_models:
times['Cloned Neural Network'] = fittime
save_times(times)
if 'Cloned Neural Network' in times:
print("Laufzeit:",times['Cloned Neural Network'],"s")
del model_clone
/usr/local/lib/python3.10/dist-packages/keras/src/engine/training.py:3103: UserWarning: You are saving your model as an HDF5 file via `model.save()`. This file format is considered legacy. We recommend using instead the native Keras format, e.g. `model.save('my_model.keras')`.
saving_api.save_model(
Laufzeit: 254.04891324043274 s
MSE_scaling_factor = (data['median_house_value'].max() - data['median_house_value'].min()) **2
scaled_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['mean_squared_error'] ]
scaled_val_loss = [math.sqrt(val * MSE_scaling_factor) for val in history.history['val_mean_squared_error'] ]
plt.figure(figsize=(20, 10))
plt.plot(scaled_loss,linestyle='--',color='lightblue')
plt.plot(scaled_val_loss,linestyle='-',color='darkblue')
plt.title('RMSE nach Epoche')
plt.xlabel('Epoche')
plt.ylabel('RMSE')
plt.grid(True)
plt.ylim(40000,70000)
plt.legend(['Training', 'Validierung'],
loc='upper left')
<matplotlib.legend.Legend at 0x7b4106ef9fc0>
feature_name='ocean_proximity'
ocean_prox_embedding = cloned_trained_model.layers[2]
ocean_prox_input = cloned_trained_model.input[0]
embedding_model_ocean_clone = keras.models.Model(inputs=ocean_prox_input, outputs=ocean_prox_embedding.output)
feature_name='county_name'
county_name_embedding = cloned_trained_model.layers[3]
county_name_input = cloned_trained_model.input[1]
embedding_model_county_clone = keras.models.Model(inputs=county_name_input, outputs=county_name_embedding.output)
Gegenüberstellung der Embedding-Grafiken:
_, axs = plt.subplots(2, 2, figsize=(20, 15))
feature_name = 'ocean_proximity'
series = x_train[feature_name]
x_coords, y_coords, labels = embedding_to_coords(embedding_model_ocean, series)
plot_embed_axis(axs[0,0], x_coords,y_coords, labels, 'Ozeannähe Embedding')
x_coords, y_coords, labels = embedding_to_coords(embedding_model_ocean_clone, series)
plot_embed_axis(axs[0,1], x_coords,y_coords, labels, 'Ozeannähe Embedding auf geklontem Model')
feature_name = 'county_name'
series = x_train[feature_name]
x_coords, y_coords, labels = embedding_to_coords(embedding_model_county, series)
plot_embed_axis(axs[1,0], x_coords,y_coords, labels, 'Landkreis Embedding')
x_coords, y_coords, labels = embedding_to_coords(embedding_model_county_clone, series)
plot_embed_axis(axs[1,1], x_coords,y_coords, labels, 'Landkreis Embedding auf geklontem Model')
plt.show()
1/1 [==============================] - 0s 28ms/step 1/1 [==============================] - 0s 45ms/step 2/2 [==============================] - 0s 5ms/step 2/2 [==============================] - 0s 6ms/step
Die Embeddings sind durch das erneute Training auf dem geklonten Modell anders, dies ist durch die randomisierten Prozesse im Training zu verantworten. Insbesondere ist es zu erwarten, dass ein gegebenes Embedding bei erneuter Durchführung gespiegelt oder gedreht sein kann, da die Relation der Punkte zueinander für das Modell maßgeblich ist, weniger die absoluten Gewichtswerte.
Dies ist bei den Embeddings für die Ozeannähe insbesondere zu beobachten. Bei den Landkreisen ist dies nicht eindeutig der Fall, vermutlich gibt es beim Training des neuronalen Netzes noch zu viele Freiheitsgrade. Es deutet sich aber an, dass gewisse Ausprägungen hier in beiden Varianten einen Extrempunkt bilden, also eher am Rand der Grafik erscheinen. Dies lässt zumindest deuten, dass die Dimensionen des Embeddings einen ähnlichen Wert repräsentieren, auf dem die entsprechenden Kategorien entsprechend stark ausschlagen.
del learning_rate
del fittime
del scaled_loss
del scaled_val_loss
del starttime
del x_train
del x_val
del y_train
del y_val
del ocean_prox_embedding
del ocean_prox_input
del county_name_embedding
del county_name_input
del feature_name
del history
del embedding_model_county, embedding_model_county_clone, embedding_model_ocean, embedding_model_ocean_clone
del series
del x_coords, y_coords
del axs
Aufgabe B-3 e)¶
Der RMSE der Validierungsdaten des Netzes aus Teilaufgabe d) ist samt Fitting-Laufzeit zu ermitteln und dem Modellvergleich hinzuzufügen. Die entsprechenden Ergebnisse der beiden Neuronalen Netz mit Embeddings sind miteinander und mit dem Neuronalen Netz aus Aufgabe B-2 b) sowie dem einfachen linearen Modell aus Aufgabe A-7 a) zu vergleichen und zu bewerten.
validation_data=sdata[data['sample']=='B'].sample(frac=1.0, random_state=seed)
x_val=validation_data.drop([target_column,'sample'],axis=1)
y_val=validation_data[target_column]
del validation_data
MSE = cloned_trained_model.evaluate([x_val['ocean_proximity'].cat.codes,x_val['county_name'].cat.codes, x_val[numeric_features].values.astype(float)], y_val.values.astype(float),verbose=0)[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmses['Cloned Neural Network'] = RMSE
print("RMSE auf Validierungsdaten: ",rmses['Cloned Neural Network'])
RMSE auf Validierungsdaten: 56769.304093019855
print("RMSEs:")
display(rmses)
print("Laufzeiten:")
display(times)
RMSEs:
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855,
'Complete Neural Network': 56852.606799089896,
'Embedding Neural Network': 56651.739797472845,
'Cloned Neural Network': 56769.304093019855}
Laufzeiten:
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956,
'Complete Neural Network': 150.87731003761292,
'Embedding Neural Network': 203.7618260383606,
'Cloned Neural Network': 254.04891324043274}
Als geklontes Netzwerk würde man eigentlich eine identische Performance als bei dem vorherigen Embedded-Modell erwarten. Dies ist allerdings nicht gegeben. Die Laufzeit weicht ebenfalls von der vorherigen Laufzeit ab. Grund dafür ist die große Rolle welche zufällige Effekte, beispielsweise die Startgewichte, im Prozess des Trainings bei neuronalen Netzwerken spielen. Dadurch kann auch die Konvergenz früher oder später einsetzen.
Es bleibt dabei, dass die Embedding-Modelle gegenüber dem Neuronalen Netz mit Encoding ("Compelte Neural Network") eine ähnliche Performance zeigen, aber höhere Laufzeiten benötigen.
Sie haben definitiv eine höhere Vorhersagekraft als das lineare Modell und sind demgegenüber fast ausschließlich zu bevorzugen, zumindest wenn es für das Training keine Echtzeitanforderungen gibt.
del MSE, RMSE
del optimizer
#für Aufgabe B-5
test_data=sdata[data['sample']=='C'].sample(frac=1.0, random_state=seed)
x_test=test_data.drop([target_column,'sample'],axis=1)
y_test=test_data[target_column]
del test_data
MSE = cloned_trained_model.evaluate([x_test['ocean_proximity'].cat.codes,x_test['county_name'].cat.codes, x_test[numeric_features].values.astype(float)], y_test.values.astype(float),verbose=0)[1]
RMSE = math.sqrt(MSE_scaling_factor*MSE)
rmsestest['Cloned Neural Network'] = RMSE
del MSE, RMSE
del MSE_scaling_factor
del esdata, sdata
del cloned_trained_model
Aufgabe B-4: Hyperparameter-Tuning von lightGBM und XGBoost¶
Aufgabe B-4 a)¶
Es ist ein lightGBM-Modell inkl. Hypertuning mit mehrdimensionalem Grid Search und vierfacher Kreuzvalidierung zu erstellen. Dabei sind die Parameter learning_rate und num_leaves zu optimieren, für den Parameter n_estimators ist der entsprechende Standardparameter des CatBoost Modells zu verwenden. Die mit dem Grid Search ermittelten besten fünf Parametersätze sind anzuzeigen und zu diskutieren. Im Anschluss ist ein Modell mit den besten Parametern auf allen Folds zu fitten, der Fehler auf den Validierungsdaten zu berechnen und die Ergebnisse dem Modellvergleich hinzuzufügen.
Vorbereiten der Daten: Die numerischen Werte müssen nicht skaliert werden, aber wir machen ein One-Hot-Encoding der Kategorien.
edata=pd.get_dummies(data, columns=categorical_features, drop_first=True)
x_train =edata[edata['sample']=='A'].drop(['sample',target_column],axis=1)
y_train=edata[edata['sample']=='A'][target_column]
x_val =edata[edata['sample']=='B'].drop(['sample',target_column],axis=1)
y_val=edata[edata['sample']=='B'][target_column]
#für Aufgabe B-5
x_test =edata[edata['sample']=='C'].drop(['sample',target_column],axis=1)
y_test=edata[edata['sample']=='C'][target_column]
Hyperparametertuning über learning_rate und num_leaves
parameter = {
'learning_rate': [0.1,0.11,0.12,0.13,0.14,0.15],
'num_leaves': [50,60,70,80,90,100]
}
starttime = time.time()
model = lgb.LGBMRegressor(force_row_wise=True,random_state=seed,verbosity=-1)
kfold = sklearn.model_selection.KFold(n_splits=4, shuffle=True, random_state=seed)
grid_search = sklearn.model_selection.GridSearchCV(estimator=model, param_grid=parameter, cv=kfold, scoring='neg_mean_squared_error', verbose=0)
grid_search.fit(x_train, y_train)
param_results = list(zip(grid_search.cv_results_['mean_test_score'], grid_search.cv_results_['params']))
param_results.sort(reverse=True)
best5=param_results[:5]
best_model = lgb.LGBMRegressor(learning_rate=grid_search.best_params_['learning_rate'], num_leaves=grid_search.best_params_['num_leaves'],force_row_wise=True)
best_model.fit(x_train, y_train)
y_pred = best_model.predict(x_val)
fittime = time.time() - starttime
times['lightGBM'] = fittime
save_times(times)
[LightGBM] [Warning] Found whitespace in feature_names, replace with underlines [LightGBM] [Info] Total Bins 1428 [LightGBM] [Info] Number of data points in the train set: 14421, number of used features: 56 [LightGBM] [Info] Start training from score 206956.332640
print("Beste 5 Parameter mit negativem MSE: ")
display(best5)
print("Beste weiter verwendete Parameter")
display(best5[0])
Beste 5 Parameter mit negativem MSE:
[(-2254056480.857634, {'learning_rate': 0.11, 'num_leaves': 70}),
(-2258871900.292492, {'learning_rate': 0.13, 'num_leaves': 50}),
(-2261010890.956952, {'learning_rate': 0.1, 'num_leaves': 60}),
(-2262200461.58615, {'learning_rate': 0.13, 'num_leaves': 60}),
(-2264451175.4730983, {'learning_rate': 0.12, 'num_leaves': 70})]
Beste weiter verwendete Parameter
(-2254056480.857634, {'learning_rate': 0.11, 'num_leaves': 70})
Die 5 besten Parameterpaare reichen in der Lernrate von 0.1 bis 0.13 und bei der Anzahl der Blätter von 50 bis 70. Insgesamt wurden Lernraten von 0.1 bis 0.15 und Blätterzahlen von 50 bis 100 getestet. Dass die Lernrate von 0.11 und Blätterzahl 70 das beste Ergebnis erzielt zeigt, dass es keine triviale Verbesserung darstellt, den Suchraum nach oben oder unten zu erweitern. Besseres Hypterparametertuning würde man wahrscheinlich eher durch Verfeinerung des Suchgitters erhalten.
rmses['lightGBM']=math.sqrt(keras.metrics.mean_squared_error(y_pred, y_val))
print("RMSE für LightGBM")
display(rmses['lightGBM'])
RMSE für LightGBM
47584.96967530819
print("Laufzeit:",times['lightGBM'],"s")
Laufzeit: 51.38854146003723 s
iterations = len(grid_search.cv_results_)
#für Aufgabe B-5:
y_pred = best_model.predict(x_test)
rmsestest['lightGBM']=math.sqrt(keras.metrics.mean_squared_error(y_pred, y_test))
del model
del kfold
del grid_search
del param_results
del best5
del best_model
del y_pred
del edata
del starttime
del fittime
del parameter
Aufgabe B-4 b)¶
Analog ist ein XGBoost-Modell inkl. Hypertuning mit mehrdimensionalem Random Search und vierfacher Kreuzvalidierung zu erstellen. Dabei sind die Parameter depth, learning_rate, subsample und colsample_bytree zu optimieren. Dabei ist die gleiche Anzahl Iterationen wie bei CatBoost und lightGBM zu verwenden. Die mit dem Randomized Search ermittelten besten fünf Parameterkombinationen sind anzuzeigen und zu diskutieren. Im Anschluss ist ein Modell mit den besten Parametern auf allen Folds zu fitten, der Fehler auf den Validierungsdaten zu berechnen und die Ergebnisse dem Modellvergleich hinzuzufügen.
Wir nehmen wieder die unskalierten aber One-Hot-Codierten Daten aus der vorherigen Aufgabe
parameter = {
'max_depth': range(3, 16),
'learning_rate': [0.02,0.04,0.06,0.08,0.1,0.12,0.14,0.16,0.18,0.2],
'subsample': [0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0],
'colsample_bytree': [0.5,0.6,0.7,0.8,0.9,1.0]
}
display(parameter)
{'max_depth': range(3, 16),
'learning_rate': [0.02, 0.04, 0.06, 0.08, 0.1, 0.12, 0.14, 0.16, 0.18, 0.2],
'subsample': [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
'colsample_bytree': [0.5, 0.6, 0.7, 0.8, 0.9, 1.0]}
Selbe Anzahl Iterationen wie bei lightGBM
print(iterations)
14
starttime = time.time()
model = xgb.XGBRegressor(random_state=seed)
kfold = sklearn.model_selection.KFold(n_splits=4, shuffle=True, random_state=seed)
random_search = sklearn.model_selection.RandomizedSearchCV(estimator=model, param_distributions=parameter, n_iter=iterations, cv=kfold, scoring='neg_mean_squared_error', random_state=seed)
random_search.fit(x_train, y_train)
param_results = list(zip(random_search.cv_results_['mean_test_score'], random_search.cv_results_['params']))
param_results.sort(reverse=True)
best5=param_results[:5]
best_model = xgb.XGBRegressor(**random_search.best_params_)
best_model.fit(x_train, y_train)
y_pred = best_model.predict(x_val)
fittime = time.time() - starttime
times['XGBoost'] = fittime
save_times(times)
print("Beste 5 Parameterkombinationen mit negativem MSE: ")
display(best5)
print("Beste weiter verwendete Parameterkombination")
display(best5[0])
Beste 5 Parameterkombinationen mit negativem MSE:
[(-2299911921.8627796,
{'subsample': 0.9,
'max_depth': 12,
'learning_rate': 0.14,
'colsample_bytree': 0.5}),
(-2307266024.92968,
{'subsample': 0.5,
'max_depth': 10,
'learning_rate': 0.06,
'colsample_bytree': 0.9}),
(-2310193416.0413957,
{'subsample': 0.7,
'max_depth': 15,
'learning_rate': 0.08,
'colsample_bytree': 0.8}),
(-2321955486.9425,
{'subsample': 0.7,
'max_depth': 13,
'learning_rate': 0.12,
'colsample_bytree': 0.9}),
(-2377892832.5986333,
{'subsample': 0.6,
'max_depth': 14,
'learning_rate': 0.16,
'colsample_bytree': 0.6})]
Beste weiter verwendete Parameterkombination
(-2299911921.8627796,
{'subsample': 0.9,
'max_depth': 12,
'learning_rate': 0.14,
'colsample_bytree': 0.5})
Die besten 5 Parametergruppen reichen bei den Teilstichproben zwischen 0.5 und 0.9, Maximaltiefe zwischen 10 und 15, Lernrate zwischen 0.06 und 0.14 und die baumweise Spaltenwahl zwischen 0.5 und 0.9. Die beste Parametergruppe liegt bis auf die Spaltenwahl innerhalb der Grenzen und auch bei der Spaltenwahl lässt sich in den Top 5 keine Monotonie erkennen, womit das Tuning auch hier nicht trivialerweise verbessert werden kann, indem die Grenzen erweitert werden.
rmses['XGBoost']=math.sqrt(sklearn.metrics.mean_squared_error(y_pred, y_val))
print("RMSE für XGBoost")
display(rmses['XGBoost'])
RMSE für XGBoost
50459.16902941104
print("Laufzeit:",times['XGBoost'],"s")
Laufzeit: 153.54116678237915 s
#für Aufgabe B-5:
y_pred = best_model.predict(x_test)
rmsestest['XGBoost']=math.sqrt(keras.metrics.mean_squared_error(y_pred, y_test))
del model
del kfold
del param_results
del best5
del best_model
del y_pred
del starttime, fittime
del iterations
del parameter
del random_search
del x_train, x_val
del y_train, y_val
del x_test, y_test
del labels
Aufgabe B-4 c)¶
Abschließend ist der Aufwand und der Nutzen der hier optimierten Modelle gegenüber dem CatBoost-Modell aus Aufgabe A-7 b) zu bewerten.
rmses
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855,
'Complete Neural Network': 56852.606799089896,
'Embedding Neural Network': 56651.739797472845,
'Cloned Neural Network': 56769.304093019855,
'lightGBM': 47584.96967530819,
'XGBoost': 50459.16902941104}
times
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956,
'Complete Neural Network': 150.87731003761292,
'Embedding Neural Network': 203.7618260383606,
'Cloned Neural Network': 254.04891324043274,
'lightGBM': 51.38854146003723,
'XGBoost': 153.54116678237915}
lightGBM hat in diesem Fall sowohl eine kürzere Laufzeit als auch ein besseres Ergebnis beim RMSE erzielt. Dies kann am gewählten Hyperparameter-Suchbereich liegen. In jedem Fall performenen beide sehr ähnlich. Von der Laufzeit her scheint lightGBM dem Namen entsprechend eine kürzere Laufzeit als XGBoost zu verzeichnen. Im Kontext dieser Aufgabe ist also in beiden Aspekten lightGBM gegenüber XGBoost zu bevorzugen.
Im Vergleich zu CatBoost hat zwar lightGBM ein leicht besseres Ergebnisse erzielt, lightGBM und XGboost benötigten dafür aber ein Hyperparameter-Tuning, während Catboost diese Leistung bereits mit Default-Parametern erzielte.
Zusammenfassend erscheint bei der Möglichkeit zum aufwendigen Tuning bei dieser Aufgabe lightGBM als bestes Modell hervorzugehen. Wenn es darum geht ohne weitere Konfiguration gute und schnelle Ergebnisse zu erzielen ist allerdings CatBoost der Favorit.
Aufgabe B-5: Finale Modellbewertung an Testdaten¶
Aufgabe B-5 a)¶
Das CatBoost-Modell aus Aufgabe A-7 sowie alle im Teil B des Notebooks enthaltenen Modelle sind abschließend mit den Testdaten (sample = C) zu bewerten und die Ergebnisse in den Testdaten abzulegen.
Die Ergebnisse auf den Testdaten wurden bereits während der vorherigen Aufgabenteile zusammengestellt:
print('RMSEs auf den Testdaten:')
rmsestest
RMSEs auf den Testdaten:
{'Complete Catboost': 45427.075729708304,
'Scaled Encoded Catboost': 46154.39938978077,
'Complete Neural Network': 55160.98627873368,
'Embedding Neural Network': 54221.61607277586,
'Cloned Neural Network': 54265.635567510915,
'lightGBM': 45852.00825481911,
'XGBoost': 48551.297191733196}
Aufgabe B-5 b):¶
Die Modellgüte der Modelle ist grafisch aufzubereiten. Welche Modelle haben die beste Prognosegüte? Gibt es im Hinblick auf Prognosegüte, Optimierungsaufwand und Laufzeit ein besonders geeignetes und empfehlenswertes Modell?
In den folgenden Grafiken bezeichnet
- "Numerical Linear" das lineare Modell auf den numerischen Testdaten
- "Categorical Linear" das lineare Modell auf den kategoriellen Testdaten
- "Complete Linear" das lineare Modell auf allen Testdaten
- "Numerical Catboost" das Catboost-Modell auf den numerischen Testdaten
- "Categorical Catboost" das Catboost-Modell auf den kategoriellen Testdaten
- "Complete Catboost" das Catboost-Modell auf allen Testdaten
- "Complete Neural Network" das neuronale Netzwerk mit One-Hot-Codierung auf allen Testdaten
- "Embedding Neural Network" das erste neuronale Netzwerk mit Embedding auf allen Testdaten
- "Cloned Neural Network" das zweite neuronale Netzwerk mit Embedding auf allen Testdaten
sorted_rmses = sorted(rmses.items(), key=lambda x: x[1])
modelle = [item[0] for item in sorted_rmses]
rmse_values = [item[1] for item in sorted_rmses]
plt.figure(figsize=(16, 8))
plt.bar(modelle, rmse_values)
plt.xlabel('Prognosemodell')
plt.ylabel('RMSE')
plt.title('RMSES der Modelle auf den Valdiierungsdaten')
plt.ylim(0.9*min(rmse_values), 1.1*max(rmse_values))
plt.xticks(rotation=45, ha='right')
plt.show()
sorted_rmses = sorted(rmsestest.items(), key=lambda x: x[1])
modelle = [item[0] for item in sorted_rmses]
rmse_values = [item[1] for item in sorted_rmses]
plt.figure(figsize=(16, 8))
plt.bar(modelle, rmse_values)
plt.xlabel('Prognosemodell')
plt.ylabel('RMSE')
plt.title('RMSES der Modelle auf den Testdaten')
plt.ylim(0.9*min(rmse_values), 1.1*max(rmse_values))
plt.xticks(rotation=45, ha='right')
plt.show()
Die auf Ensemble-Trees basierten Modelle schneiden insgesamt besser ab als die Modelle der neuronalen Netze. Es ist nicht auszuschließen, dass die neuronalen Netze durch feinere Justierung und mehr Tuning- und Trainingsaufwand besser abschneiden würden. Unter den Ensemble-Baum-Modellen heben sich Out-of-the-box Catboost und lightGBM besonders durch einen niedrigen Fehler hevor. Auch hier ist es wieder möglich, dass XGBoost nur aufgrund unseres Hyperparameter-Tunings zurückbleibt. In Angesicht der Tatsache, dass für lightGBM zunächst ein Hyperparameter-Tuning notwendig war ist hier CatBoost mit Default-Parametern das besonders empfehlenswerte Modell.
del sorted_rmses
del rmse_values
del modelle
Teil C: Clustering unter Verwendung von Feuergefährdungszonen¶
Aufgabe C-1: Erweiterung des Hauptdatensatzes um Feuergefährdungszonen¶
Aufgabe C-1 a)¶
Der Datensatz california_fire_zones.shp ist einzulesen (siehe Anhang 1), im Anschluss sind 10 zufällige Zeilen davon auszugeben. Die absolute Häufigkeitsverteilung der Werte des Merkmals FHSZ ist zu erstellen. Weiter ist eine Funktion zu schreiben, welche für ein gegebenes Polygon die Anzahl an Eckpunkten zurückgibt. Unter Verwendung dieser Funktion ist ein Histogramm zur Anzahl an Eckpunkten der Polygone des eingelesenen Datensatzes zu erstellen, bei welchem die x-Achse logarithmiert dargestellt wird, und zu kommentieren. Schließlich ist ein Plot der Feuergefährdungszonen zu erzeugen und mit einer Karte von Kalifornien zu hinterlegen.
Einlesen der Daten und Umkodierung der Ortsdarstellung via Längen- und Breitengrade
firezones = gpd.read_file(folder_path+'california_fire_zones.shp').to_crs(epsg=4326)
firezones.sample(10, random_state=seed)
| FHSZ | FHSZ_Descr | Shape_Leng | Shape_Area | geometry | |
|---|---|---|---|---|---|
| 1092 | 3 | Very High | 2.289474e+04 | 6.115388e+06 | POLYGON ((-119.53744 37.30656, -119.54267 37.3... |
| 468 | 2 | High | 2.194811e+04 | 5.640945e+06 | POLYGON ((-117.94386 36.16922, -117.94392 36.1... |
| 1332 | 3 | Very High | 1.382022e+04 | 5.107807e+06 | POLYGON ((-120.98480 40.58090, -120.98452 40.5... |
| 827 | 2 | High | 1.926041e+04 | 4.366668e+06 | POLYGON ((-121.30206 41.37308, -121.30207 41.3... |
| 929 | 3 | Very High | 1.574675e+04 | 3.785108e+06 | POLYGON ((-117.82993 34.06540, -117.82984 34.0... |
| 699 | 2 | High | 1.643052e+05 | 2.190773e+07 | POLYGON ((-121.51306 39.60514, -121.51212 39.6... |
| 377 | 2 | High | 1.192624e+05 | 1.240135e+07 | POLYGON ((-117.04513 33.48070, -117.04514 33.4... |
| 1066 | 3 | Very High | 1.923844e+06 | 2.043624e+09 | POLYGON ((-121.63036 36.58334, -121.63034 36.5... |
| 447 | 2 | High | 3.594848e+04 | 6.195179e+06 | POLYGON ((-118.47049 35.64034, -118.47050 35.6... |
| 1454 | 3 | Very High | 1.017415e+05 | 3.337294e+07 | POLYGON ((-120.54643 41.53049, -120.54679 41.5... |
firezones['FHSZ'].value_counts()
FHSZ 3 627 2 519 1 361 Name: count, dtype: int64
Anzahl der Ecken per Polygon zusammenstellen
cornercounts = firezones['geometry'].apply(ncornerpoints)
cornercounts.describe()
count 1507.000000 mean 2018.330458 std 6613.219157 min 7.000000 25% 267.500000 50% 618.000000 75% 1295.500000 max 107983.000000 Name: geometry, dtype: float64
plt.figure(figsize=(15, 10))
plt.hist(cornercounts, bins='auto')
plt.xscale('log')
plt.xlabel('Anzahl an Eckpunkten')
plt.ylabel('Häufigkeit')
plt.title('Histogram der Eckpunkte (automatische Topfbreite gemäß Freedman-Diaconis)')
plt.show()
logbins = np.logspace(np.log10(cornercounts.min()), np.log10(cornercounts.max()), num=15)
plt.figure(figsize=(15, 10))
plt.hist(cornercounts, bins=logbins)
plt.xscale('log')
plt.xlabel('Anzahl an Eckpunkten')
plt.ylabel('Häufigkeit')
plt.title('Histogram der Eckpunkte (exponentiell wachsende Töpfe)')
plt.show()
Die Anzahl der Ecken der Polygone schwank stark, von 7 bis über 100000. Die Mehrzahl der Polygone hat weniger als 1000 Ecken.
Mit steigender Anzahl an Eckpunkten gibt es immer weniger Polygone, welche diese repräsentieren. Wenn man die Töpfe exponentiell wachsen lässt zeigt sich eine Art Normalvertielung mit Peak bei ca. 500, d.h. der Logarithmus der Polygonecken ist annähernd normalverteilt.
del cornercounts
del logbins
california_img_path = folder_path+'california.png'
cal_img = plt.imread(california_img_path)
legend_values = [1, 2, 3]
firezones.plot(column='FHSZ', figsize=(20, 10),legend=True,cmap=plt.cm.hot_r,vmin=-0.5, vmax=4, alpha=0.8,legend_kwds={'ticks': legend_values, 'label': 'Feuergefährdungszone'})
plt.xlabel('Längengrad')
plt.ylabel('Breitengrad')
plt.title('Feuergefährdungszonen')
axisrange = [-124.55, -113.95, 32.45, 42.05]
plt.axis(axisrange)
plt.imshow(cal_img, extent=axisrange)
plt.show()
del axisrange, legend_values
Aufgabe C-1 b):¶
Der Hauptdatensatz, der nach der Bearbeitung von Teil A vorliegt, ist um die in Teilaufgabe a) eingelesenen Feuergefährdungszonen zu erweitern. Zu diesem Zweck soll festgestellt werden, ob und ggf. in welcher Feuergefährdungszone die Longitude-Latitude-Paare der einzelnen Zeilen liegen. Das Resultat ist in der neuen Spalte fire_zone abzuspeichern. Liegt ein Longitude-Latitude-Paar in keiner Feuergefährdungszone, so ist als Wert 0 einzutragen. Liegt ein Longitude-Latitude-Paar in mehr als einer Feuergefährdungszone, so ist eine aktuariell angemessene Lösungsstrategie zu diesem Problem durchzuführen. Schließlich soll das Ergebnis überprüft werden, indem ein geeigneter Scatter-Plot der Daten erzeugt wird.
data_geometry = gpd.GeoDataFrame(geometry=gpd.points_from_xy( data['longitude'], data['latitude']),crs=firezones.crs)
data_geometry['data_index']= data.index
matchingfirezones = gpd.sjoin(data_geometry, firezones, how="left", predicate="within")
matchingfirezones.reset_index(drop=True)
| geometry | data_index | index_right | FHSZ | FHSZ_Descr | Shape_Leng | Shape_Area | |
|---|---|---|---|---|---|---|---|
| 0 | POINT (-117.04000 32.54000) | 0 | NaN | NaN | NaN | NaN | NaN |
| 1 | POINT (-117.09000 32.55000) | 1 | NaN | NaN | NaN | NaN | NaN |
| 2 | POINT (-117.06000 32.55000) | 2 | NaN | NaN | NaN | NaN | NaN |
| 3 | POINT (-117.04000 32.55000) | 3 | NaN | NaN | NaN | NaN | NaN |
| 4 | POINT (-117.12000 32.56000) | 4 | NaN | NaN | NaN | NaN | NaN |
| ... | ... | ... | ... | ... | ... | ... | ... |
| 20667 | POINT (-121.93000 41.86000) | 20635 | NaN | NaN | NaN | NaN | NaN |
| 20668 | POINT (-123.83000 41.88000) | 20636 | NaN | NaN | NaN | NaN | NaN |
| 20669 | POINT (-124.16000 41.92000) | 20637 | NaN | NaN | NaN | NaN | NaN |
| 20670 | POINT (-124.14000 41.95000) | 20638 | 1506.0 | 3.0 | Very High | 8.494362e+04 | 6.690529e+07 |
| 20671 | POINT (-122.64000 41.95000) | 20639 | 1503.0 | 3.0 | Very High | 1.385404e+06 | 9.765650e+08 |
20672 rows × 7 columns
print("Anzahl Datenpunkte:",len(data_geometry))
print("Feurzonen-Überlappungen:",len(matchingfirezones))
Anzahl Datenpunkte: 20640 Feurzonen-Überlappungen: 20672
Es gibt also Datenpunkte, welche in mehreren Feuerzonen liegen.
Bei mehreren Zonen pro Datenpunkt wählen wir die höchste Gefährdungsstufe
data_to_firezone = matchingfirezones.groupby('data_index')['FHSZ'].max().fillna(0).rename('fire_zone', inplace=True)
fdata = pd.merge(data, data_to_firezone, left_index=True, right_on='data_index')
print("Anzahl Datenpunkte mit Feuerzone:",len(fdata))
Anzahl Datenpunkte mit Feuerzone: 20640
Scatterplot aller Datensätze, mit farblicher Codierung der Feuergefährdung
# Datensatz sortieren, damit der Plot immer die höchste Feuerzone zeigt
fdata_sorted=fdata.sort_values(by='fire_zone')
plot_axis =fdata_sorted.plot(
kind="scatter", x="longitude", y="latitude",
s=fdata_sorted["population"] / 100,
xlabel="Longitude", ylabel="Latitude",
label="Population",
c=fdata_sorted['fire_zone'], cmap=plt.cm.hot_r,vmin=-0.5, vmax=4,
colorbar=True,
legend=True,
sharex=False, figsize=(20, 10))
plot_axis.collections[0].colorbar.set_label('Feuergefährdungszone')
legend_values = [0, 1, 2, 3]
plot_axis.collections[0].colorbar.set_ticks(legend_values)
axisrange = [-124.55, -113.95, 32.45, 42.05]
plt.axis(axisrange)
plt.imshow(cal_img, extent=axisrange)
plt.show()
del fdata_sorted, legend_values
Man erkennt hier eindeutig die gleiche Farbverteilung wie im Plot der Feuergefährdungszonen
del data_geometry, data_to_firezone
del firezones
del matchingfirezones
del plot_axis
Aufgabe C-2: Einsatz von Clustering bei der Tarifierung von Feuerversicherung¶
Aufgabe C-2 a)¶
Um den aus Aufgabe C-1 erhaltenen Datensatz zu clustern, ist der k-Means-Algorithmus anzuwenden. In das Clustering sollen lediglich die Merkmale latitude, longitude, median_house_value und fire_zone einfließen. Die Wahl der Cluster-Anzahl ist (etwa durch eine geeignete Heuristik) zu begründen. Ein Scatter-Plot zur Visualisierung der Cluster ist zu erzeugen.
cdata = fdata[['latitude', 'longitude','median_house_value','fire_zone']]
Wir plotten die Silhouettenkoeffizienten um eine Idee für eine gute Clusteranzahl zu erhalten
candidates = range(2, 11)
scores = []
for candidate in candidates:
kmeans = sklearn.cluster.KMeans(n_clusters=candidate, random_state=seed, n_init='auto')
cluster_labels = kmeans.fit_predict(cdata)
score = sklearn.metrics.silhouette_score(cdata, cluster_labels)
scores.append(score)
plt.plot(candidates, scores, marker='o')
plt.xlabel('Clusteranzahl')
plt.ylabel('Silhouettenkoeffizient')
plt.title('Silhouettenkoeffizienten')
plt.show()
del score, scores
2 und 5 scheinen gute Kandidaten für die Anzahl der Cluster zu sein. Um mehrere Aspekte der Daten abbilden zu können wählen wir 5 Cluster.
numclusters = 5
kmeans = sklearn.cluster.KMeans(n_clusters=numclusters, n_init='auto',random_state=seed)
kmeans.fit(cdata)
kmeans_labels = kmeans.labels_
Die Datenpunkte werden auf einer Karte mit ihrer Clusterzugehörigkeit farblich dargestellt. Die Größe der Punkte stellt den jeweiligen Median-Hauswert dar.
plot_axis =fdata.plot(
kind="scatter", x="longitude", y="latitude",
s=fdata["median_house_value"] / 10000,
xlabel="Longitude", ylabel="Latitude",
label="Median-Hauswert",
c=kmeans_labels, cmap='rainbow',
colorbar=True,
legend=True,
sharex=False, figsize=(20, 10))
plot_axis.collections[0].colorbar.set_label('Cluster')
axisrange = [-124.55, -113.95, 32.45, 42.05]
plt.axis(axisrange)
plt.imshow(cal_img, extent=axisrange)
plt.show()
del candidate, candidates
del cluster_labels
del kmeans, kmeans_labels
Aufgabe C-2 b)¶
Um den aus Aufgabe C-1 erhaltenen Datensatz zu clustern, ist das Verfahren Agglomerative Hierarchical Clustering anzuwenden. In das Clustering sollen lediglich die Merkmale latitude, longitude, median_house_value und fire_zone einfließen. Die Wahl der Cluster-Anzahl ist (etwa durch eine geeignete Heuristik) zu begründen. Ein Scatter-Plot zur Visualisierung der Cluster ist zu erzeugen.
Wir nutzen wie im vorherigen Aufgabenteil wieder 5 Cluster, basierend auf dem Silhouettenkoeffizienten
np.random.seed(seed)
agglomerative = sklearn.cluster.AgglomerativeClustering(n_clusters=numclusters)
agglo_labels = agglomerative.fit_predict(cdata)
plot_axis =fdata.plot(
kind="scatter", x="longitude", y="latitude",
s=fdata["median_house_value"] / 10000,
xlabel="Longitude", ylabel="Latitude",
label="Median-Hauswert",
c=agglo_labels, cmap='rainbow',
colorbar=True,
legend=True,
sharex=False, figsize=(20, 10))
plot_axis.collections[0].colorbar.set_label('Cluster')
axisrange = [-124.55, -113.95, 32.45, 42.05]
plt.axis(axisrange)
plt.imshow(cal_img, extent=axisrange)
plt.show()
del agglo_labels
del agglomerative
del axisrange
del cal_img
del california_img_path
del fdata, cdata
del numclusters
del plot_axis
Aufgabe C-2 c)¶
Die Clusterings der beiden vorherigen Teilaufgaben sind zu vergleichen, indem auf drei Gemeinsamkeiten oder Unterschiede eingegangen wird.
Übereinstimmendes Clustering im Inland
Die beiden Clusterverfahren erzeugen in der breiten Fläche des Inlands (welches sich durch vergleichsweise niedrige Median-Hauswerte auszeichnet) fast identische Clusterzuordnungen (Farben Türkis und Orange). Im vergleich mit den Feuergefährdungszonen scheint das orangene Cluster tendenziell eine hohe Feuergefährdung und Türkis eine niedrige zu haben.
Vergleichbare Kategorien in den urbanen Küstengebieten
Das Clustering in den urbanen Küstegebieten (mittlere bis hohe Hauspreise) führt zumindest zu vergleichbaren Clustern. Eine grobe (vereinfachte) Interpetation der übrigen 3 Cluster wären
- A: mittlere Preise, erhöhte Feuergefahr
- B: hohe Preise, niedrige Feuergefahr
- C: hohe Preise, hohe Feuergefahr
Diese sind in den Plots dargestellt mit
- A: KMeans lila, Agglo grün
- B: KMeans grün, Agglo rot
- C: KMeans rot, Agglo lila
Unterschiedliche Einteilung im Detail
Im Detail liefern beide Algorithmen jedoch an einigen Stellen eine unterschiedliche Zuteilung zu den Clustern. Ein Beispiel wären die beiden Insel-Datenpunkte bei ca. (-118.5, 33.5) Die beiden Algorithmen haben in ihren Ergebnissen also vergleichbare Kategorien gewählt, aber unterschiedliche Kriterien und Grenzen für die konkrete Zuordnung.
Aufgabe C-2 d)¶
Ein Ansatz, wie die von einem der Clustering-Verfahren erzeugten Cluster bei der Tarifierung von Feuerversicherung in Kalifornien eingesetzt werden können, ist in zwei bis drei Sätzen zu erläutern.
Das Clustering kann eingesetzt werden, um für eine Feuerversicherung eine geringe Anzahl an unterschiedlichen Preiskategorien festzulegen, die sich anschließend pro Landkreis zuordnen lassen. Die höheren Versicherungsleistungen bei hochwertigeren Häusern des Gebiets und die Wahrscheinlichkeit eines Brandschadens mittels der Feuergefährdungszone können dabei zentrale Tarifmerkmale sein. Das Miteinbeziehen des geografischen Standpunktes kann hier einen "Smoothing"-Effekt bewirken, sodass die Zuordnungen zu den Preisstufen in einem Gebiet relativ einheitlich sind.
Teil D: GPU-Ausführung und Vergleiche¶
Aufgabe D-1: Ergebnisse des CPU-Notebooks persistieren¶
Aufgabe D-1 a)¶
display(times)
{'Complete Linear': 0.11070752143859863,
'Complete CatBoost': 12.42766809463501,
'Scaled Encoded CatBoost': 6.519219160079956,
'Complete Neural Network': 150.87731003761292,
'Embedding Neural Network': 203.7618260383606,
'Cloned Neural Network': 254.04891324043274,
'lightGBM': 51.38854146003723,
'XGBoost': 153.54116678237915}
if retrain_models:
timesdf=pd.DataFrame(times.items(), columns=['Modellname','Laufzeit'])
timescsvpath=results_path+'times.csv'
timesdf.to_csv(timescsvpath, index=False)
del timesdf, timescsvpath
display(rmses)
{'Numeric Linear': 74508.93978488311,
'Numeric Catboost': 48912.79736396997,
'Categorical Linear': 91193.2953907611,
'Categorical Catboost': 89735.68699792186,
'Complete Linear': 69098.04473488631,
'Complete Catboost': 48361.638797917374,
'Scaled Encoded Catboost': 48406.630856312855,
'Complete Neural Network': 56852.606799089896,
'Embedding Neural Network': 56651.739797472845,
'Cloned Neural Network': 56769.304093019855,
'lightGBM': 47584.96967530819,
'XGBoost': 50459.16902941104}
rmsesdf=pd.DataFrame(rmses.items(), columns=['Modellname','RMSE auf Validierungsdaten'])
rmsescsvpath=results_path+'rmses.csv'
rmsesdf.to_csv(rmsescsvpath, index=False)
del rmsesdf, rmsescsvpath
display(rmsestest)
{'Complete Catboost': 45427.075729708304,
'Scaled Encoded Catboost': 46154.39938978077,
'Complete Neural Network': 55160.98627873368,
'Embedding Neural Network': 54221.61607277586,
'Cloned Neural Network': 54265.635567510915,
'lightGBM': 45852.00825481911,
'XGBoost': 48551.297191733196}
rmsestestdf=pd.DataFrame(rmsestest.items(), columns=['Modellname','RMSE auf Testdaten'])
rmsestestcsvpath=results_path+'rmsestest.csv'
rmsestestdf.to_csv(rmsestestcsvpath, index=False)
del rmsestestdf, rmsestestcsvpath
CPU-Hardware:
psutil.cpu_count()
2
!cat /proc/cpuinfo
processor : 0 vendor_id : GenuineIntel cpu family : 6 model : 79 model name : Intel(R) Xeon(R) CPU @ 2.20GHz stepping : 0 microcode : 0xffffffff cpu MHz : 2199.998 cache size : 56320 KB physical id : 0 siblings : 2 core id : 0 cpu cores : 1 apicid : 0 initial apicid : 0 fpu : yes fpu_exception : yes cpuid level : 13 wp : yes flags : fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush mmx fxsr sse sse2 ss ht syscall nx pdpe1gb rdtscp lm constant_tsc rep_good nopl xtopology nonstop_tsc cpuid tsc_known_freq pni pclmulqdq ssse3 fma cx16 pcid sse4_1 sse4_2 x2apic movbe popcnt aes xsave avx f16c rdrand hypervisor lahf_lm abm 3dnowprefetch invpcid_single ssbd ibrs ibpb stibp fsgsbase tsc_adjust bmi1 hle avx2 smep bmi2 erms invpcid rtm rdseed adx smap xsaveopt arat md_clear arch_capabilities bugs : cpu_meltdown spectre_v1 spectre_v2 spec_store_bypass l1tf mds swapgs taa mmio_stale_data retbleed bogomips : 4399.99 clflush size : 64 cache_alignment : 64 address sizes : 46 bits physical, 48 bits virtual power management: processor : 1 vendor_id : GenuineIntel cpu family : 6 model : 79 model name : Intel(R) Xeon(R) CPU @ 2.20GHz stepping : 0 microcode : 0xffffffff cpu MHz : 2199.998 cache size : 56320 KB physical id : 0 siblings : 2 core id : 0 cpu cores : 1 apicid : 1 initial apicid : 1 fpu : yes fpu_exception : yes cpuid level : 13 wp : yes flags : fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush mmx fxsr sse sse2 ss ht syscall nx pdpe1gb rdtscp lm constant_tsc rep_good nopl xtopology nonstop_tsc cpuid tsc_known_freq pni pclmulqdq ssse3 fma cx16 pcid sse4_1 sse4_2 x2apic movbe popcnt aes xsave avx f16c rdrand hypervisor lahf_lm abm 3dnowprefetch invpcid_single ssbd ibrs ibpb stibp fsgsbase tsc_adjust bmi1 hle avx2 smep bmi2 erms invpcid rtm rdseed adx smap xsaveopt arat md_clear arch_capabilities bugs : cpu_meltdown spectre_v1 spectre_v2 spec_store_bypass l1tf mds swapgs taa mmio_stale_data retbleed bogomips : 4399.99 clflush size : 64 cache_alignment : 64 address sizes : 46 bits physical, 48 bits virtual power management:
psutil.virtual_memory()
svmem(total=13609451520, available=11391655936, percent=16.3, used=1930776576, free=9702719488, active=809914368, inactive=2841772032, buffers=250802176, cached=1725153280, shared=11325440, slab=151568384)